Pa x, ¢ . 13)) .. . if. t)! 59?.»49 ' O '§u£?{ RF I: z :1 ch! r: 44““ I..3D)\.l .A . .rtlrr . . “Md”... 5:: 4 fr}.1r;lf: .‘vt‘lvx. a vireo‘lefi) PIJDYOei t p A . .Iv‘v 9 l.)-..v.. fifiuxu , . 7 v .e. n 45.. .... .. . . . ‘ . ‘ . . V . . I w) ‘5 CF. Ilu.r.. .. : . . , it.(.u 5.23‘ 15:5 . . f . 1 . ‘ if t. Havana»... Z. . .l . ar Sr ( tics-r .fr Kraut...” (LI ”’1, . . u I .. ‘ , .:.t.1tvyu.v , L .. . , t4! 4: 13!.) OF ‘1’! L. . . } . .l .E v Isl- le-n ybtl..r .gnvlo livl (. . , 43.53. ‘1 ‘ .. tfpvtozi (.Ei>.ll)4 by .14.. ..Jltt.... rky u. L..t.\i:. 35" It 51: .1! an (.71.. 2 ' x 5 . : . . .. ., . ‘ . . . H.,...r..\ .s; l x x . .. . . . A! $.23.» .) Faun...‘ . . . I o 2 . I513. . 011.1 4 Mus. C.‘ in‘ . . .. C; )3): r»!!! ."lr It“... 0’ if. ) . . : vttnrvlalto! 0!; ab‘rv a; z. 119????) LI: 1-. ‘ ..:v.oam)...(t.l))!\.4$?r.lp ’01:!» (19:10:); t.£ln!...\$; f .5 .‘fvlfliltl,’ 11,51... .uciri 1 I . . . In... ..£.\v€cr.r x- . . v .1 . . ‘3.) ‘ . ...‘.‘ 1!}, i 1.13 fifty. .5. . ., ‘ v . . , re,“ '~rA€§u . .. ‘ 3.-.". sll. u if; . .. 1 A . . o .. . \ .I . I... .yl 4 (I. g. . rior 4»... .. t. L J - . fl . . u..| . , 5.. . .1." 1 5 I r, . x5... rl, :1111251. tr .3 Jr "(thaw-PL.) I: . . . . . » .vlli‘ . . , I ‘uf. (xv ...|...flov “M. L 3‘3.-. Z... . 1:-Llc£ .,..‘ 1» III, ' 1a .. 1.: ct .. . yfi... \ , (I. ll '1’ fizhirt . ‘ f4. 1. a: F . 5' .I' _ .. . r! )‘u ‘DH ‘1'“?! (‘ wjbfl7945 MICHI IGAN STA ..;___ IIIHIIIIIIII/IIIIIIIHIIIMU [3111771 L LIBRARYO Michigan Sign LUnivordty This is to certify that the dissertation entitled Modelling and Control of Mechanical Systems with Stick-Slip Friction presented by Stephen Charles Southward has been accepted towards fulfillment of the requirements for Ph . D . Mechanical Engineering degree in 4 / MS U is an Affirmative Action/Equal Opportunity Institution 0-12771 PLACE IN RETURN BOX to remove thie checkout from your record. TO AVOID FINES return on or before date due. DATE DUE DATE DUE DATE DUE ll Kb . 1'“ [[1113] L l” .3 l 4."; ‘3 -auhf APR 06 1993 , ’1 24 .T .r l 2 [—7 MSU Is An Affirmative Adieu/Equal Opportunity Inetitmion Modelling and Control of Mechanical Systems With Stick-Slip Friction By Stephen Charles Southward A DISSERTATION Submitted to Michigan State University in partial fulfillment of the requirements for the degree of DOCTOR OF PHILOSOPHY Department of Mechanical Engineering 1990 1;") Ln C) Q) Abstract Modelling and Control of Mechanical Systems with Stick-Slip Friction By Stephen Charles Southward Whenever two mechanical surfaces in a dynamical system are forced together in direct contact, a friction force is induced at the point of contact which acts to resist the relative motion of the surfaces. The result is often a stick-slip friction force that is highly nonlinear, and generally has an undesirable influence on the system dynamics. Its primary function is to remove useful energy from the system. With a representative model of this phenomena, control schemes can be developed and evaluated to compensate for the undesirable effects. In this dissertation, a model is presented which, when incorporated into a dynamic system model, will provide the stick-slip dynamics observed in practice. Using this model, classical control techniques are evaluated for their effectiveness. A Proportional+Derivative (PD) control law is the simplest classical control law to regulate the position of a one-degree-of-freedom (l-DOF) system. A stable set of multiple equilibrium points is found to exist for the 1-DOF system under PD control, and the steady-state error is guaranteed to be bounded. Integral control (PID) is added in an attempt to remove the steady-state error, but certain types of slipping force models are shown to promote the generation of limit cycles. Neither of these classical techniques are able to effectively regulate the position of the 1-DOF system. A nonlinear compensation force for stick-slip friction is developed to supplement a PD control law applied to the 1-DOF mechanical system. The choice of a discontinuous compensation force is motivated by the requirement that the desired reference be a unique equilibrium point of the system. The stick-slip friction force, modelled with a sticking force and a slipping force, generates discontinuous state derivatives. A Lyapunov function is introduced to prove global asymptotic stability of the desired reference using a modification of the direct method for discontinuous systems. Stability is verified numerically as well as experimentally. .The nonlinear compensation force is robust with respect to the character of the slipping force which is assumed to lie within piecewise linear bounds. Exact knowledge of the static friction force levels is not required, only upper bounds for these levels. Copyright COPyright by Stephen Charles Southward 1990 Certification This is to certify that the dissertation entitled Modelling and Control of Mechanical Systems with Stick-Slip Friction presented by / ,/ , . W, I mm m J Am] 331960 8 1 p *n Charles Southward . I Date DOCTOR OF PHILOSOPHY Candidate has been accepted towards fulfillment of the requirements for a DOCTOR OF PHILOSOPHY DEGREE in Mechanical Engineering at Michigan State University ; ; Date M r‘ \VKMHQSL’VNA AQV‘I\ 3094“) Charles R. MacCluer Date Professor, Dept. of Mathematics Z/Zaw / 9ko I770 0 R0 enberg Date Professor, Dept. of Mech Engineering 3947—3; 93,996) MW» Fathi Salam Date Associate ofessor, Dept. of Electrical Engineering flIfi/V, /3C’ IDaZVtJ Charles Seebeck Professor, Dept. of Mathematics U 5%.; 14m? 22 {270 7 7 Date St‘e en W. Shaw Associate Professor, Dept. of Mechanical Engineering Dedication To my parents, Charles and Mildred Southward. vi Acknowledgements I would like to extend my most sincere gratitude to my academic advisor Dr. Clark J. Radcliffe. He has taught me much more than I had anticipated learning through this experience. Special thanks go to Dr. Chuck MacCluer for his advice and assistance through our brief association. I doubt that I could have finished this dissertation without his guidance and expertise. Both of these men shall remain a source of positive inspiration to me in my future pursuits, academic or otherwise. I am very grateful to claim them as friends and mentors. Each member of my guidance committee played an integral part in the development of this research effort, and I would like to extend my gratitude to each of them. To Dr. Steven Shaw for his support and guidance with the geometrical analysis. To Dr. Ronald Rosenberg for his assistance with interpretation of the nonlinear simulation results. And finally to Dr. Fathi Salam and Dr. Charles Seebeck for providing me with the fundamental expertise in each of their respective disciplines. I am proud to have been able to work with such an outstanding group of individuals. Special thanks are extended to Bob Rose for his mechanical design expertise. I would also like to thank the General Motors Corporation for providing funds to equip the Dynamic Systems Laboratory, and the CASE Center for CAEM for providing computational services. Finally, I would like to thank all my friends at Michigan State, especially Andy (Cap’n DPH Albano) Hull, who has truly shared this experience with me. vii Table of Contents List of Tables List of Figures Chapter 1 - Introduction 1.1 The Study of Stick-Slip Friction 1.2 The Control of Systems with Stick-Slip Friction 1.3 Chapter Overview Chapter 2 - Mathematical System Models 2.1 The Physics of Stick-Slip Friction 2.2 The Stick-Slip Friction Model 2.3 The Reaction Force 2.4 Common Slipping Force Models 2.5 The One-Degree-Of-Freedom System Model Chapter 3 - PD Control Response 3.1 Theoretical Investigation 3.2 Experimental Investigation Chapter 4 - PID Control Response 4.1 Theoretical Investigation 4.2 Simulated Time Responses 4.3 Limit Cycle Simulations 4.4 Experimental Investigation Chapter 5 - Robust Nonlinear Compensation 5.1 Derivation of the Control Law viii \OQQUIUJ 12 15 19 21 21 23 29 29 34 41 47 52 52 5.2 Surface of Discontinuity and State Trajectories 55 5.3 Stability Theorem and Proof 57 5.4 Simulation Results 61 5.5 Experimental Results 68 Chapter 6 - Conclusions 72 6.1 Summary 72 6.2 Recommendations for future work 75 Appendix A - Simulation Models 76 A.l Numerical Simulation 76 A2 Stick-Slip Friction Model 76 A3 One-Degree-Of-Freedom Model 78 Appendix B - Experimental Controller 82 B.1 Controller Implementation 82 B.2 PID Controller 82 B.3 PD Control with Nonlinear Stick-Slip Compensation 94 8.4 Common Control Modules 106 References l 10 ix Table 3.1 3.2 4.1 4.2 4.3 4.4 4.5 4.6 4.7 4.8 4.9 4.10 4.11 5.1 5.2 5.3 5.4- 5.5 List of Tables Static friction bound results from PD control experiment. Experimental parameters for PD control. Common parameters for all three PID simulations. PID Simulation parameters for the Coulomb model. PID Simulation parameters for the stiction model. PID Simulation parameters for the exponential model. Classifications for each of the three stick-slip models. Common parameters for PID limit cycle simulations. Numerical results of fitted iteration functions for the Coulomb model. Numerical results of fitted iteration functions for the stiction model. Numerical results of fitted iteration functions for the exponential model. Experimental parameters for PID control. Slipping force parameters fitted to experimental limit cycle response. Common parameters for all nonlinear compensation simulations. Simulation parameters for nonlinear compensation with the Coulomb model. Simulation parameters for nonlinear compensation with the stiction model. Simulation parameters for nonlinear compensation with the exponential model. Experimental parameters for nonlinear compensation. SA Page 26 27 36 37 38 42 43 45 46 48 50 62 63 65 69 List of Figures Figure 2.1 The sticking force model with unequal static friction levels. 2.2 Example nonlinear slipping force model bounded by linear bounds. 2. 3 Coulomb plus viscous friction model. 2.4 Stiction plus viscous friction model. 2.5 Exponential plus viscous friction model. 2.6 The 1-DOF conceptual stick-slip mass system. 3.1 Typical PD phase portrait and the set of multiple equilibria E PD. 3 .2 Experimental 1-DOF system schematic diagram. 3.3 Experimental brake to alter the friction characteristics. 3 .4 Experimental system time response for PD control. 4.1 Parallel trajectories inside the sticking band S. 4.2 Typical PID trajectory with sticking on both sides of the origin. 4. 3 PID Control simulated time response with the Coulomb model. 4.4 PID Control simulated time response with the stiction model. 4.5 PID Control simulated time response with the exponential model. 4. 6 Iteration functions for PID control with the Coulomb model. 4.7 Iteration functions for PID control with the stiction model. 4. 8 Iteration functions for PID control with the exponential model. 4.9 Experimental system time response for PID control. 4.10 Simulated limit cycle response fitted to experiment. 5.1 The nonlinear proportional feedback control force. xi Page 11 11 16 17 18 19 22 24 27 32 34 37 39 40 45 46 49 5 1 5.2 5.3 5.4 5.5 5.6 5.7 Nonlinear addendum, g(x), to the Lyapunov function candidate. Simulated time response for nonlinear compensation with Coulomb model. Simulated time response for nonlinear compensation with stiction model. Simulated time response for nonlinear compensation with exponential model. Experimental time response for nonlinear compensation with braking. Experimental time response for nonlinear compensation without braking. xii 58 .65 67 7O 71 Chapter 1 Introduction 1.1 The Study of Stick-Slip Friction Ever since man began using machines to improve his control over nature, stick-slip friction has been present to hinder and sometimes even help in this process. Stick-slip friction, sometimes known as dry friction, is a natural resistance to relative motion between two contacting bodies. This resistance can be reduced to effect the prevention of surface wear with the use of lubricants, as the Egyptians, Greeks, and Romans were well aware (Davison, 1957). Although Aristotle recognized and wrote of the reality of friction during his time, it was nearly 2000 years before a quantitative study of friction took place (Bowden and Tabor, 1964). In the mid-fifteenth century, Leonardo da Vinci enunciated the first definitive laws of friction based on direct experiments illustrated in his notebooks. The concept of a coefficient of friction is clearly expressed in his statement, “The friction of a polished smooth surface resists the engine with a force equal to one quarter of its weight” (V inci, circa 1500c). Although the details were not quite correct, the fundamental concepts behind his laws of friction still hold me today. It is remarkable to note that these laws were developed 200 years before Newton presented a clear definition of force. The rediscovery of these laws by the French engineer Guillaume Amontons was presented in 1699 to the Royal Academy of Sciences in France, but was met with some skepticism. Amontons, working with dirty surfaces, determined the coefficient of friction to be one third for all surfaces. These observations were verified by Charles Augustine Coulomb in 1781 working with clean surfaces, and his results were also presented to the Royal Academy. Through meticulous experimental work, Coulomb observed that static friction is usually higher than kinetic friction, but for metals, the static friction is almost identical to the kinetic friction. With strong encouragement from the Royal Academy, Coulomb pioneered the research effort to determine the fundamental causes of friction. Much of the friction research through the nineteenth century was concerned with the mode of energy dissipation in friction. It is now accepted that friction is largely due to interfacial adhesion and to the energy expended in deforming the surfaces (Bowden and Tabor, 1964). Although this research continues to some degree today, friction research has branched off into many new disciplines. The areas of primary interest to this research are modelling stick-slip friction, and the associated effects of these models on the control of dynamic systems. Many models have been formulated for this type of friction, but only a few are regularly used in practice today. The analytical simplicity of these common models are based on the fundamental research of Coulomb. Through experimental evidence, Coulomb observed that kinetic friction was nearly independent of the relative velocity of the smfaces, thus the Coulomb friction model was postulated. For some surfaces, the kinetic friction was found to be lower than the static levels. This model is currently accepted as the stiction model. A linear viscous model, where the friction force is proportional to the velocity, emerged from investigations with various lubricants. Linear combinations of these models are usually able to predict observed system dynamics, and because of their simplicity, they are the most common models found in practice. Extensive research efforts have been directed toward characterizing the actual friction force in experimental systems. Rabinowicz (1951) investigated the nature of the transition region, from static to kinetic friction, when stationary metal surfaces were set into morion. He postulated that the transition might be a continuous drop in contrast to the discontinuous transition of the stiction model. Walrath (1984) introduced a first-order dynamic friction model which effectively produced a continuous transition from static to kinetic friction. In order to more accurately represent the stick-slip dynamics observed in practice, Karnopp (1985) introduced a model with explicit simulation considerations. Kolston (198 8) utilized electrical circuit analysis equivalences to study a mechanical system with stick-slip friction. Research efforts have also been directed toward identifying the associated model parameters for an experimental system. Based on the in-phase and quadrature power dissipated when exciting a normal mode, Tomlinson and Hibbert (1979) developed a technique to evaluate the magnitude of the nonlinear friction force and the hysteretic damping constant. Parameter identification has been posed as an optimization problem, and Cheok et al. (1988) employed a modified simplex algorithm to identify the system parameters for servomechanisms with stick-slip friction. Motivated by the desire to make our machines more efficient and/or more effective, the direct effect of stick-slip friction on dynamical systems has also been studied extensively. Den Hartog (1930) presented the exact solution to the forced response of a one-degreeof-freedom (l-DOF) system using a linear combination of coulomb and viscous friction. Shaw (1986) extended this work to include a stiction model using bifurcation theory to study the nonlinear dynamics. Gogoussis and Donath (1987) examined the fundamental issues associated with stick-slip friction, and their consequences on the design and performance of robots. 1.2 The Control of Systems with Stick-Slip Friction Feedback control techniques implemented for precise positioning of mechanical systems is very difficult because stick-slip friction is discontinuous when the relative velocity is zero. In order to use the well-developed design tools from linear control theory, dynamic system models are typically linearized about some operating point. Even with an accurate model, attempts to linearize or neglect stick-slip friction at zero velocity are either impossible or produce invalid predictions of system response. For this reason, regulation of position is much more difficult than tracking a desired trajectory. This observation is implicit in the literature where the regulator problem is investigated almost exclusively. The effects of nonlinear stick-slip friction on feedback control were examined using a vector graphic technique (T ustin, 1947), which is similar to the modern describing function analysis technique. Tou and Schultheiss (1953) thoroughly studied the impact on control with static and kinetic friction using the describing function technique. More recently, Townsend and Salisbury (1987) investigated the effect of coulomb friction and stiction on force control through a compliant transmission with integral feedback. The degree of difficulty of the control problem is clearly indicated with the types of problems that are attacked in the literature. The regulation of l-DOF systems with stick- slip friction is studied almost exclusively. This may be due in part to the large number of l-DOF regulator applications requiring stick-slip friction compensation. It may also reflect the problems associated with the control of multi-DOF systems which have stick-slip friction. A broad range of compensation techniques have been developed for an even broader range of l-DOF applications. Gilbart and Winston (1974) developed a model reference adaptive control system to compensate for bearing friction in an optical tracking telescope. Walrath (1984) used a new dynamic friction model as the basis for a predictive, adaptive friction compensator, and demonstrated a significant improvement in the performance of airborne servo mechanisms. Kubo et a1. (1986) presented a feedback compensation controller for a robot manipulator which attempts to eliminate the effect of the static friction force. Canudas' et al. (1987) developed an adaptive control scheme to feedback linearize a DC motor system with nonlinear friction. Yang and Tomizuka (1987) developed and verified an adaptive pulse-width control scheme which adaptively updates the width and height of the applied input control force pulse to regulate each of the X-Y table coordinates. Armstrong (1988) utilized experimental results of friction behavior to precompute motion torques for open loop compensation of a brush type DC servo motor with gearing. In this dissertation, we develop asimple nonlinear compensation technique for the regulation of a 1-DOF system with stick-slip fiiction. The choice for the compensation force is motivated by the requirement that the desired reference be a unique equilibrium point of the system. This leads to a discontinuous nonlinear proportional feedback force which appears to be a bang-bang force in a region near the desired reference. This compensation force is robust to the character of the slipping friction force which is assumed to lie within piecewise linear bounds. Only upper bounds on the static friction levels are required to obtain stability with this control law. A modified form of Lyapunov’s direct method is used to prove global asymptotic stability of the equilibrium point. Stability is illustrated with a numerical simulation, and experimental test results verify the stability of the reference position. 1.3 Chapter Overview The physical phenomena of stick-slip motion is examined in Chapter 2. A general mathematical model is presented which produces the desired stick-slip dynamics when incorporated into a dynamic system model. This model is compared with the classical friction models found in practice today. The state equations for the one-degree—of-freedom system are also presented in this chapter. Chapter 3 investigates the theoretical time response of the l-DOF system with stick- slip friction under PD control. A procedure for estimating the static friction levels is derived from the theoretical analysis of this system. Experimental static friction levels are determined with this procedure, and the theoretical time response is verified with an experimental system. The response of the 1-DOF system with stick-slip fiiction under PID control is explored in Chapter 4. A theoretical analysis of the PID time response indicates that trajectories may have a limit cycle behavior. Further numerical analysis verifies the existence of limit cycle behavior for some slipping force models. An experimental time response verifies the predicted theoretical limit cycle response. A robust nonlinear compensation force is derived in Chapter 5. This compensation supplements a standard PD control law to regulate the position of the 1-DOF system with stick-slip friction. Global asymptotic stability is proven with Lyapunov stability methods. Simulation and experimental time responses verify the global stability of the origin. Chapter 6 concludes and summarizes this dissertation with recommendations for future work. Chapter 2 Mathematical System Models it is necessary to know the nature of the contact which this weight has with the smooth surface where it produces friction by its movement, because different bodies have difi‘erent kinds of friction; because if there shall be two bodies with difi’erent surfaces, that is that one is soft and polished and well greased or soaped, and it is moved upon a smooth surface of a similar kind, it will move much more easily than that which has been made rough by the use of lime or a rasping-file. ' Leonardo da Vinci Forster Manuscript S.K.M. 11, 86v. and 87r. A model is required to represent the physical characteristics of stick-slip friction. In this chapter, the qualitative character of this type of fiiction is investigated, and a general mathematical model is presented. This mathematical stick-slip friction model, when incorporated into a dynamical system model, will exhibit the same behavior as that observed in the corresponding physical system. Common qualitative slipping force models, found in practice today, and their mathematical representations are reviewed here. The dynamical system equations for a one-degree-of-freedom system with stick-slip friction are also presented in this chapter. 2.1 The Physics of Stick-Slip Friction When two mechanical surfaces are in contact with each other, a stick-slip friction force is generated at the contacting surface. This stick-slip friction force has two defining characteristics which can effectively be separated into two force components. A slipping, or kinetic friction force dissipates energy when the relative velocity of the surfaces is non- zero. The slipping force thus impedes the relative motion of the surfaces. A sticking, or static friction force acts when the relative velocity is zero, therefore it does not dissipate energy. The sticking force has the ability to keep the relative velocity zero under certain limiting conditions. These two forces do not act simultaneously though they are dependent on each other as will be shown. The static friction force will take on whatever value is necessary to keep the relative velocity zero, limited by positive and negative extremes. These are the smallest forces necessary to start a relative motion in either direction. Leonardo da Vinci first discovered that these maximum static friction forces generally followed two empirical laws. The first law; “Friction produces double the amount of effort if the weight be doubled” (V inci, circa 1500a), expresses the proportionality of the static friction force to the normal force between the surfaces (Halliday and Resnick, 1978). The second law is; “The friction made by the same weight will be of equal resistance at the beginning of the movement although the contact may be of different breadths or lengths” (V inci, circa 1500b). This is a statement of the observation that the static friction forces are approximately independent of the surface area of contact. The static coefficient of friction is defined to be the ratio of the maximum static friction force to the normal force. It is difficult in practice to obtain a globally constant coefficient of friction due to the strong dependence on local surface properties (Bowden and Tabor, 1954). The slipping force is generally a function of the relative velocity between the surfaces when that velocity is non-zero. As the static friction force bounds are dependent on the normal force between the surfaces, so too is the slipping force. In practice, a kinetic coefficient of friction is determined as a function of velocity. The slipping force is the product of the normal force and the kinetic coefficient of friction. For many dry polished surfaces of dissimilar contacting materials, the kinetic coefficient of friction is approximately constant (Amontons, 1699; Euler, 1750; Coulomb, 1785). Experimental determination of the kinetic coefficient of friction is even more difficult than the static coefficient (Bowden and Tabor, 1954). The analysis in this dissertation will be restricted to systems where the normal force between the surfaces is constant. With this simplification, the static friction force bounds are constant, and the slipping force can be represented by a function of velocity only. Even with this simplification, the stick-slip friction force remains highly nonlinear. The overall stick-slip friction force can be expressed as a function of relative velocity when that velocity is non-zero, but must be determined by an algebraic equilibrium constraint when zero. The apparent switching of input-output causality on the stick-slip friction force causes a real modelling problem. One set of dynamic equations are required for the slipping motion, and another set for the sticking motion. A model for representing the stick-slip phenomena with constant causality was inu'oduced by Karnopp (1985). With this method, only one set of state equations are required. 2.2 The Stick-Slip Friction Model The model for the stick-slip friction force 1",, used throughout this analysis is a modification of the model presented by Karnopp (1985). This model allows for the evaluation of the friction force during sticking and slipping motions. It also eliminates the . numerical problems associated with near-zero-velocity force computation. True zero velocity will rarely be numerically attained, therefore an artificial zero region is defined around zero velocity [—Ot,0t]. The conceptual stick-slip friction model is achieved with equation (2.1) only in the limit as a—>0. F, = Fw(v)[7t(v)] + mexl - am} (2.1a) x _ 1 lvl>a >0 21b M" o lvlsa ’a ’ (‘ ) where v is the relative velocity between the surfaces, and F is a reaction force which depends on the particular system configuration. The reaction force is, briefly, that force which is required to keep the relative velocity zero and is the solution of an algebraic 10 constraint which depends on the overall system equations. A thorough description of this reaction force can be found in the following section, and in Karnopp (1985). Throughout the following analytic investigation, the artificial zero parameter a is taken to be identically zero. Non-zero values of or are only used in simulations to insure that the numerical integration algorithms remain stable. Any velocity which is numerically within a of zero is taken to be zero velocity. The sticking function, Fm» provides values of the friction force at zero velocity. This term is used to determine whether the contacting surfaces will stick, or break free from the static friction forces. The positive and negative limits on the static friction forces are given by ll",+ and F: respectively. As mentioned above, the static friction force bounds are assumed to be constant throughout this analysis, but they are not presumed to be equal in magnitude. F," F 2 F: > 0 Fm,(F) = F F: < F < F: (2.2) If F S F: < 0 When the relative velocity is zero, v=0, and the reaction force F is between the static friction levels, the damping mechanism provides a force which keeps the surfaces from slipping. The surfaces remain stuck until F is greater in magnitude than the respective static friction force (Figure 2.1). Figure 2.1. The sticking force model with unequal static friction levels. The slipping function, F . ,hp, provides values of the friction force at non-zero velocity, and is represented by two functions of the relative velocity, line) = E*(v-a)u(v)+B‘(v+a)u(-v). (2.3) where u(.) is the right-continuous Heavyside step function. The function F; (.) defines the slipping force for positive velocity only, and the function F:(.) defines the slipping force for negative velocity only. Two separate functions are required for this representation since the slipping force is not assumed symmetric. An example slipping force with (1:0 is shown in Figure 2.2. Figure 2.2. Example nonlinear slipping force model bounded by linear bounds. 12 Throughout the following analysis, F; (.) and F: (.) are assumed to be arbitrary continuous functions which satisfy a Lipschitz condition (Vidyasagar, 1978, pp. 79-80) over their respective domains. The slipping force is assumed to dissipate energy at all non-zero velocities, and therefore is bounded within the first and third quadrants as indicated in Figure 2.2. We assume there exist constants b1 2 b0 > 0, F: S F: < 0, and F? 2 F”+ > 0, which define piecewise linear bounds for the slipping functions (Figure 2.2). bov S Fflv) S F;++ blv Vv > 0 (2.4a) If + blv S F;(v) S bov Vv < 0 (2.4b) As the velocity approaches zero from the right or left, the limiting values of the slipping force do not, in general, exceed the respective static friction force levels, P? and If (Rabinowicz, 1951), thus we have natural lower limits on the magnitudes of F? and 15;”. This property will be required in the following analysis, but it should not be taken as a restriction since this property is commonly observed in practice. 2.3 The Reaction Force The reaction force F is the sticking force which is required to keep the relative velocity zero, limited by the static friction bounds. This force is the solution of an algebraic equilibrium relation which depends on the particular system equations. Presented here is an algorithm for computing the reaction force(s) for a general N-DOF mechanical system. The overall mechanical system can contain rotational and/or uanslational degrees of freedom. This restriction is made since these types of motions have relative surface motions and thus stick-slip friction exists at the surface contact. A general notation is required to set up the framework for this algorithm. For the N-DOF system, there are N generalized coordinates given by the vector q , N associated 13 generalized velocities given by the vector 6;, and N associated stick-slip friction forces given by the vector Fd . There are also M generalized inputs given by the vector u. q = lap-".407 (2.5a) it drip-".4517 (2.5b) F, = [Fd,,---,Fm]r (2.5c) u = [%,...,uM]T (25d) Note that only N stick-slip friction forces are allowed in this system representation. This should not be taken as a restriction since multiple stick-slip friction forces associated with a single degree of freedom can always be combined, or lumped together, to form a single stick-slip friction force. The individual components of the stick-slip friction vector are expressed similar to the model in equation (2.1). Since this is a multi-degree of freedom representation, the relative velocity between two surfaces associated with the generalized coordinates is expressed as the difference between two generalized velocities. F “’(v(‘))[7t(v“’)] + F§L(F)[l- -l(v“’)] (2.6a) where the relative velocities are 1);": (qj— q‘,). j,k e {0,--- ,N} (2.6b) The relative velocity associated with the i‘“ generalized coordinate is the difference between the j‘” and k‘h generalized velocities. Since it is possible for the slipping velocity to be relative to ground, or zero velocity, the 0 index is used to represent this situation (i.e. q'0 '=‘ 0). The particular indices (j,k) associated with the i‘h generalized coordinate are determined from the system model. 14 Associated with each component of the stick-slip friction force then is a slipping force which is a function of a particular relative velocity, and a sticking force which .is a function of its associated reaction force. These reaction forces only need to be computed when the associated relative velocity becomes zero. With the following general state equations, ii = f(q.<'r.F,.u) (2.7) the algorithm for the computation of the reaction forces can explicitly be stated. The first step is to compute the N relative velocities associated with each of the stick-slip friction forces. For every relative velocity that is non-zero, evaluation of the associated stick-slip friction force is explicit. The slipping force can be evaluated at this velocity to generate the pr0per friction force. v5.2 at 0 => F4, = 17,2265?) (2.8) When the i‘’1 relative velocity is zero, the sticking force will take on values to keep that velocity zero, limited by the associated static friction bounds for that generalized coordinate. Initially, we will assume that the sticking force is not limited, therefore from equation (2.2) v5.2 = 0 => 1;. = F; (2.9) In order to keep the relative velocity zero, there must be no relative acceleration, thus (61' -ék)l = 0 => fj(qvq’ Fd’u) = ft(Q9q9Fd:u) (210) The resulting algebraic equation (2.10) has the unknown reaction force F} appearing explicitly in the stick-slip friction force vector. For every relative velocity which is zero, there will exist an associated algebraic equation. There will always be an equivalent 15 number of equations and unknowns which can be solved simultaneously for the reaction forces. Repeating the above process, the unlimited reaction forces are obtained through the simultaneous solution of the algebraic equilibrium constraints, and finally they must be limited by their respective sticking force functions. Fr = 17.22413) (2.11) Under the premise that the sticking forces are not limited, when a particular relative velocity becomes zero, it must remain zero for all time thereafter. The physical lirrritations on the sticking forces indicate that the relative velocity may not remain zero. 2.4 Common Slipping Force Models Many slipping force models are available for use in dynamic simulations or analysis (Kolston, 1988; Craig, 1988; Canudas, et al., 1987; Gogoussis and Donath, 1987; Shaw, 1986). Most of the standard slipping force models found in the literature can be classified into three categories. In each of these categories, a linear viscous damping term is included since it is the most common damping model, and the simplest model for analysis. The first category is a Coulomb plus viscous damping model. The slipping force approaches the respective static friction bound in a linear fashion from the left or right (Figure 2.3). l6 Figure 2.3. Coulomb plus viscous friction model. The slipping force equations for this model are F;(v) = F," + b‘v (2.12a) F;(v) = F; + b’v, (2.12b) where b*>0 and b'>0 are the viscous damping coefficients for positive and negative velocity respectively. The slipping force can thus have different viscous damping characteristics for positive and negative velocity. The next category is a stiction plus viscous damping model. The slipping force approaches a kinetic friction level which is at a lower level than the respective static friction bound (Figure 2.4). For many materials the kinetic friction level is approximately half of the static friction level (Miller, 1977). Figure 2.4. Stiction plus viscous friction model. The slipping force equations for this model are F;+(v) = F? + b+v (2.13a) F;(v) = F: + b'v, (2.13b) where b*>0 and b'>0 are the viscous damping coefficients for positive and negative velocity respectively, and F,“ > F: > 0, and F: < F: < 0. The last category is an exponential model with viscous damping. This slipping friction force exponentially decays from the static level to a kinetic friction force with increasing velocity (Tustin, 1947; Armstrong, 1988). After a certain velocity, the viscous term dominates and the slipping force increases linearly with velocity (Figure 2.5). Figure 2.5. Exponential plus viscous friction model. The slipping force equations for this model are F;(v) = Fj+(at-1«;*)[1—e“"'3’]+b+v (2.14a) F;(v) = F; + (F; — 17;)[1 - e“"'”] + b‘v, (2.14b) where b*>0 and b’>0 are the viscous damping coefficients for positive and negative velocity respectively, and Ii" > if > 0, and F,“ < F: < 0. The two velocity constants, v: and v; , are the characteristic velocities where 63% of the drop from the static to the kinetic force level occurs. Though the limiting values as velocity approaches zero in (2.14) are the static friction bounds, a more general exponential model can be obtained by lowering these limits below the static levels. Any pair of contacting surfaces will have a slipping force relation which may be approximately represented by one of the three models presented above. In a purely analytic investigation, the primary motivation for the choice of slipping force model is simplicity. If an experimental system is to be modelled, a complex model may be required for a more accurate response prediction. It is also reasonable to expect that some physical mechanical systems will exhibit a slipping force characteristic which can be represented by the simplest slipping force model. 19 2.5 The One-Degree-Of-Freedom System Model The l—DOF system under investigation is a mass constrained to move in one dimension along a horizontal ground surface. Displacement and velocity of the mass are measured relative to the ground. Stick-slip friction is present between the mass and the supporting surface (Figure 2.6). Displacement, Velocity x,v Friction Force Mas Control Force a \! => . Figure 2.6. The l-DOF conceptual stick-slip mass system. The two first-order state equations for this system are given by x = v (2.15a) v = ()4,)(F-F;,). (2.15b) These equations can also be used to represent a rotational l-DOF system with the proper change of units. Linear displacements and velocities must be replaced with rotational displacements and velocities, forces must be replaced with torques, and the mass must be replaced with a rotational inertia. However, if the stick-slip friction model presented in section 2.2 above is to be used for a rotational l-DOF system, the assumption of a constant ~ normal force between the contacting surfaces must remain valid. The differential equations of motion (2.15) are linear except for the stick—slip friction force, F, . The stick-slip friction force will be completely defined when a particular slipping force model is chosen and an expression for the reaction force is obtained from the system equations. In this example, the relative velocity is just the absolute velocity of the mass (with respect to ground). In order to keep the mass at rest when the velocity is zero, the acceleration must remain zero. Setting the acceleration term in (2.15b) to zero yields the 20 algebraic relation whereby the reaction force is clearly equal to the input control force F. As long as the applied force F remains between the static friction bounds, the mass will remain stuck. When the applied force is greater than the respective static friction level, there will be a non-zero acceleration, and the mass will begin to slip. Chapter 3 PD Control Response The very rapid friction of two thick bodies produces fire. Leonardo da Vinci Codex Atlanticus We would like to find a simple control law to regulate the position of the mass under the influence of stick-slip friction. Classical control theory offers several attractive possibilities. The simplest classical control technique is a Proportional+Derivative (PD) law, where the control force is proportional to the position error and also to its derivative. A substantial amount of useful information can be obtained through the study of this control law applied to the mass system. A bounded steady-state error is found to exist for the 1-DOF system with a PD control law. The theoretical l-DOF system response to a PD control law is analyzed and the results are experimentally verified in this chapter. 3.1 Theoretical Investigation A PD control law is simple and easy to implement. When applied to the l-DOF mass system (Figure 2.6) with a type of stick-slip friction modelled by (2.1-2.4), state trajectories end up at or near the desired reference with a bounded steady-state error (Hahn, 1967; Shaw, 1986). The steady-state error bounds can be explicitly determined by solving for the equilibrium points of the system. Without loss of generality, the origin of the phase space is taken as the desired reference position. The PD control force then simplifies to F = —K,;c -— Kdv, (3-1) where K P > 0 is the proportional gain, and K d > 0 is the derivative gain. Solving the 1- DOF system equations (2.15) with (3.1) for the equilibrium points 21 N M i=0 => v=0 (3.2) v = 0 => [—Kpr—Kdv-EEPMvPFWJl -7t.(v))] = 0 => -K,x—Fm=0. (3.3) Note that since the velocity is zero at equilibrium, Fm will take on any value of the applied force limited by the static friction levels. In the equilibrium constraint (3.3), the applied force is just -K,;r. There are two maximum displacements from the origin, xL and x”, whereby the proportional control forces are equal to the respective static friction forces. 'Ihese maximum displacements are the steady-state error bounds, and are given by F“ F ‘ xL = -616} x” = -[K—’} (3.4) From (3.2-3.4), we see that there exists a set of multiple equilibrium points which can be written as E”, = {(x,v)lv = 0,xL S x S x,,}. The existence of this set of multiple equilibrium points, E”), (Figure 3.1) is undesirable. Ax'=v f5 32 Figure 3.1. Typical PD phase portrait and the set of multiple equilibria E». The PD-space equilibrium set is an invariant set and Lyapunov's direct method for stability of invariant sets can be employed to verify the stability of Em (Hahn, 1963, and 23 1967). If the PD control gains are chosen to stabilize the undamped system, then with the proper choice of a Lyapunov function, it can be shown that the equilibrium set E”, is globally asymptotically stable. Every trajectory will terminate with zero velocity and some bounded steady-state position error. More importantly, no limit cycles can exist in the PD control space for damping functions of the form presented in Chapter 2 with the given assumptions. There are only two trajectories in the PD space that take the mass position directly to the origin. Unless initial conditions are specified on one of these two trajectories, there will be some non-zero bounded steady-state error. The steady-state error bounds (3.4) can be reduced by increasing the proportional gain K P, but at the expense of large control forces for states far from the origin. For any finite proportional gain, there will exist a non-zero steady-state error. This suggests that in order to get zero steady-state error, the proportional gain must be infinite at the origin, or the control force may need to be discontinuous. 3.2 Experimental Investigation An experimental l-DOF system was available to validate the predicted theoretical system reSponses. A cart, constrained to move in one dimension along a track, was cable driven by a DC servo-motor. A control voltage to the DC servomotor was transformed to a force which could push the cart in either direction along a one meter track. A rotational potentiometer, rigidly coupled to the drive-train, admits electrical measurement of cart positions. The general experimental system schematic is shown in Figure 3.2. Position Sensor DC Servomotor Figure 3.2. Experimental l-DOF system schematic diagram. The friction present in the l-DOF experimental system was actually distributed throughout the various components in Figure 3.2. For example, bearing friction existed in the servomotor, cart wheels, position sensor, and drive cable shafts. Because it is assumed to be a l-DOF system, the various friction components could be lumped into one global stick-slip friction force. In order to increase the number of types of friction possible for investigation, a friction brake was installed. With this brake, the overall friction force could be altered to provide a range of static friction levels, and also vary the slipping force characteristics. Tightening Screws Drive Shaft :2 i .4 Base Support .\. IIIIIII Figure 3.3. Experimental brake to alter the friction characteristics. 25 A full-scale drawing of the experimental friction brake is shown in Figure 3.3. The spring loaded brake pressed down on the drive shaft of the cable drive system. By changing the compression in the springs and using different brake materials, different desirable or undesirable friction characteristics could be obtained. All control algorithms investigated here were implemented digitally using 12-bit A/D and D/A converters to interface with the position sensor and DC motor respectively. The static calibration constants for the actuator and sensor were determined to be 0.75936 N/V, and 19.895 V/m respectively (Southward, 1986). Based on these calibration constants and the digital hardware configuration, the resolution of position measurements was 0.000245 m, and the resolution of force output was 0.0037 N. These resolutions allow ranges of approximately 10.5 m. for position, and 21:76 N. for force. The velocity of the cart, which is needed to compute the derivative control force, was estimated with a finite difference algorithm using position data. First order digital filters were available for both the position data and the estimated velocity data. The controller algorithm was implemented on a DEC LSI-11/23 computer in Fortran IV. Sampling rates up to 200 Hz. were possible with the floating-point implementation (Appendix B). Equation (3.4) suggests a method of obtaining estimates for the upper bounds on the static friction forces. Assuming the static friction bounds are constant (i.e. independent of position), the PD control law can be applied to the experimental system to determine the sticking limits xL and x” for a given proportional gain K P. Theoretically, the cart must stick at all points between the sticking limits. These limits were experimentally determined by starting the PD control law with the mass inside the sticking zone, and manually pushing it toward either bound. Beyond the sticking limit, the proportional force is greater than the static friction level, and the mass will be forced back into the sticking zone. 26 This test was performed on the experimental l-DOF system with PD control and the friction brake applied to the drive shaft. Over the indicated range of proportional gains, with the derivative gain set equal to zero, the static friction bounds were determined from the measm'ed sticking limits (Table 3.1). Table 3.1. Static friction bound results from PD control experiment. KAN/In) E” (N) F. (N) 8 1.768 -1483 10 1.791 -1479 12 ' 1.705 -1431 14 1.761 -1372 16 1.750 -1449 20 1.765 -1.442 25 1.710 -1.416 30 1.713 -1435 35 1.776 -1.468 40 1.705 -1420 As seen from the above results, the static friction bounds appear to vary slightly with position. The mean values and respective standard deviations from the above data are: Ff: 1.744 :1: 0.033 N. and Ff: -1.439 i 0.033 N. From the relatively small standard deviations, it is reasonable to assume that the static friction bounds are not strongly position dependent, thus we are justified in assuming it to be constant for the ‘braked’ system. The existence of the sticking zone can now be verified with an experimental time response. Using the estimates for the static friction forces, the sticking zone limits are estimated for a given proportional gain. These computed sticking zone bounds are shown in Table 3.2 with all the pertinent control parameters for the experimental setup. 27 Table 3.2. Experimental parameters for PD control. Mass m l 375 kL. PD Conuol K, 20.0 N/m Gains K4 0.5 Ns/m. Static Friction Ff 1.744 N Bounds P: -1.439 N. Sticking xL -0.0872 m. Limits xL’ 0.0719 m. Initial 10,) 0.3 m. Conditions V0.) 0.0 m/s. As indicated above, the mass is initially at rest outside of the sticldng zone. The time responses of the mass position and the applied control force are shown in Figure 3.4, thus demonstrating that a non-zero equilibrium can exist for PD control. 4 -4 23 x” 0: — - - - - - - - —— -2; l I x1. -4: i 5 F -6- -1 : x,xL,x,,: (x10 m.) F: (N.) -3 ....... , ....... ,...44.,....... 0.0 0.5 1.0 1.5 2.0 Figure 3.4. Experimental system time response for PD condo]. Voltages from the sensor and actuator were sampled at 2000 Hz. per channel and 1000 samples were taken from each to generate the experimental time response. From Figure 28 3.4, the mass gets inside the sticking zone but has enough momentum to carry it a bit further (see Figure 3.1). After one second the mass is stuck at a non-zero position, and a non-zero control force is still being applied. The control force is approximately 1.0 N (Figure 3.4), which is less than the estimated static friction level therefore we would expect the mass to remain stuck. Chapter 4 PID Control Response Define to me why one who slides on the ice does not fall. Leonardo da Vinci Codex Atlanticus The classical method for removing steady-state error from a system is to include integral action with the PD control. Thus a Proportional+Integral+Derivative (PID) control force is proportional to the position error, its derivative, and its integral. Though the PID control law is not the simplest, it is a logical choice for attempting to remove the undesirable steady—state error observed with the PD law. The effects of PID control on the 1-DOF system response are investigated in this chapter. A theoretical analysis demonstrates that state trajectories appear to have a limit cycle character. Further investigation with numerical simulations provides evidence that limit cycles can exist under PID control. A property of slipping force models is shown to promote limit cycle generation. This limit cycle behavior is experimentally verified. Based on the time responses of the l-DOF system under P11) control for the three classes of slipping force models presented in Chapter 2, a qualitative characterization of the slipping force present in the experimental system is determined. 4.1 Theoretical Investigation The inclusion of an integral force term to the control law introduces an added dimension to the state space. To simplify the following PID analysis, the static friction bounds are assumed to be constant and equal in magnitude, therefore F? = —F;' a F}. The stick-slip friction force is also assumed to be symmetric for analytical simplicity. The general system equation for PID control is 29 30 mt + F, = er + Kdé + KJ" edt, (4.1) where the actuating error, e = x, — 1:. Again, for the regulator problem, and without loss of generality the desired states are chosen to be at the origin (i.e. xd(t) = 0 and v4(t) = 0 Vt 2 to). This system is also linear with the exception of the stick-slip friction force F4. The sysrem equations can be written in state-space canonical form for x40) = 0 as - r P 1 y 0 1 0 ”y 0 x = 0 0 l x - O (4.2) v -(fr) {5.2) {E} v (Q) _ m m m _b t. m . with the identities v = x and y = if: edt. The stick-slip friction force is the same as that presented in Chapter 2, (2.1-2.4), except for the evaluation of the reaction force. For the PID controller, the reaction force is again equal to the control force which is now given by F=—Kdv—Krt-K‘.y. (4.3) This PID {control force is used in the evaluation of the sticking force (2.2). The PID control solution space for (4.2) is three dimensional with state coordinates (y,x,v). There are several subspaces in this phase space which have interesting features. The first subspace of interest is the equilibrium set, which is the set of all points in the space where the derivatives are zero. The solution of (4.2) with zero derivatives yields y=0 =9 x=0 (4.4) i=0 => v=0 (4.5) 9:0 => [-KJ—pr-Kdv-th(v)-F¢,(l-h(v))]=0 => -Ky-F;m=0. (4.6) 31 From the definition of the sticking force (2.2), the system is in equilibrium whenever IKyI S 15;. Two new constants, yL and y” , can be defined similar to those found in Chapter 3 for the steady-state error bounds on displacements. The set of equilibrium points for the PID phase space can then be written as: Em, = {(y,x,v)| x = 0, v = 0, yL S y S y"), where F: _ E Yr. ‘ “[E} yH - [K] (4-7) Notice the difference between this set and the set of equilibrium points for the PD control space, E”). Even though it is not a single point set, any trajectories which end up in Em, reach both the desired states x4 = 0 and v4 = 0. Unlike the stable PD equilibrium set. this equilibrium set is unstable in the sense that a small perturbation from an equilibrium point may result in convergence to a new equilibrium or divergence from the set entirely. As will be shown below, for certain slipping force models, trajectories near the set Em, diverge away from the set E ”D into a stable limit cycle. There exists a subspace where the system dynamics are reduced to first order independent of the type of slipping force. For trajectories to exist in the x-y plane, the velocity must be zero, v = 0, and there must be zero acceleration. From (4.2) and the definition of the stick-slip friction force, zero acceleration is possible whenever IKy + K pgtl S F}. Using the previous definitions for the coordinate bounds (3.4) and (4.7), one can write the following equivalent expressions {YLS (y+(KP/Ki)x)Sy,, v: 0 st (x+(K,/K,)y)sx,,' (4.8) These inequalities define a band in the x-y plane shown in Figure 4.1. This “sticking” band can be defined in a way that is similar to the definition of the equilibrium set as the union of two sets; S = S1 U52, where S1 = [(y,v,x)l yL S (y+ (Kp/ Ki)x) S y”, v = 0,): > 0}, and S2 = [(y,v,x)l yL S (y + (K P / K,)x) S )p, v = 0, x < 0}. Inside this band, the system 32 dynamics (4.2) are reduced to first order in the state y. Trajectories in S represent “sticking” trajectories where the mass remains stationary until the force due to the integral of the error builds up to a level greater than the maximum static friction force. Given an initial condition (ywxmvm) e S at time t” , state solutions are given by y,(t) x“(t — t”) + y” x,(t) = x t” S t S tr v,(t) 0 The mass will break free at time t4 given by it -(y.. + (K, / Kart...) x” > 0 t4 = [n+(1/x”) yL-(yn+(Kr/K‘)x’°) I” < 0 Trajectories in S are parallel to the y-axis, and their direction is dependent on the sign of x, (Figure 4.1). These trajectories show the equilibrium, Em, to be unstable for small perturbations in x and y. Figure 4.1. Parallel trajectories inside the sticking band S. 33 The remainder of the three-dimensional PID control space is divided up into two regions N and P, characterized by the sign of the velocity component of state trajectories inside each region. These regions are: P = P,UP,, where P, = {(y,v,x)l y, 2 (y+(Kp/K,)x), v = 0}, and P2 = {(y,v,x)lv > 0}. The other half of the space is similarly defined as: N = N, U N2 , where N, = {(y,v,x)ly,, S (y+ (K P/K,)x),v = 0} , and N2 = {(y,v,x)lv < 0}. From (2.2) and (4.3), we have: (y,x,v) e N, = (dv/dt) S 0 and (y,x,v) e P, z: (dv/dt) 2 0. All initial conditions (yp,vp,xp) e P at time t" have trajectories with non-negative velocity component until time tn, , and all initial conditions (y”,vm,x~,) e N at time t” have trajectories with non-positive velocity component until time t”. The final times are the first time after t" or t” where the velocity is zero, i.e. in P: vP(tp,) = 0 and vp(t) at 0 Vt e (twtfi), and in N : v,(t,,) = 0 and vn(t) at 0 v: e (tw,t,,). This sticking band S plays an important role in the possible generation of limit cycles in the PID space. The x-coordinate of a point in S physically represents the displacement of the mass away from the origin (the desired position) where the mass remains stuck until the integral force term breaks it free. Near the origin, all trajectories in S, enter N and then enter 5,. Also near the origin, u'ajectories in S2 enter P and similarly re- enter S,. Physically, the mass is stuck on one side of the origin until the force due to the integral term builds up enough to break free from the maximum static friction force. It then overshoots the origin and sticks on the other side. A similar process is repeated on the return motion. A typical trajectory is shown in Figure 4.2. 34 Figure 4.2. Typical PID trajectory with sticking on both sides of the origin. Since the general stick-slip friction force is nonlinear, and only piecewise linear for a few special cases, an explicit solution is not possible. It is not possible to determine analytically whether a limit cycle exists for a general slipping force function. 4.2 Simulated Time Responses Numerical simulations of the time responses of the 1-DOF system under PID control with the various slipping force models introduced in Chapter 2 provides a means of establishing the existence or non-existence of limit cycles. Using the stick-slip friction model with the three common classes of slipping force models, the simulated time responses are evaluated for comparison. All simulations were performed in double precision using the DIFFEQ program (Southward, 1989) on a Digital Equipment Corporation Micro-VAX computer (Appendix A). Several integration algorithms were utilized to verify that each would yield the same prediction. The algorithms tested were an Euler method, a fourth order fixed step size Runge-Kutta method, and a fourth order Adams variable step size predictor-corrector 35 method. The final results were produced from a Hamming predictor-corrector general (HPCG) method, with an initial step size of 0.001 s. The differential equations of motion from Chapter 2 (2.15) are expressed in continuous time and were integrated as such. However, in a digital controller implementation, the sensor measurements and the actuator output are quantized into discrete levels. In order to obtain a more realistic simulation, quantized states (simulating the measurement process) are used to compute the control force, which is in turn quantized (simulating the output process). Signals from the position and velocity sensors are assumed to be quantized by 12-bit ND and D/A converters. Based on the experimental parameters from Chapter 3, it is reasonable to assume resolutions of 0.00025 m. for positions, and 0.005 N. for the force. These resolutions will allow ranges of $0.5 m., and i100 N. It is also reasonable to assume a resolution of 0.01 m/s. for velocities, providing a range of i200 m/s. Using the above assumptions, the following parameters were chosen to approximately model an experimental system (T able 4.1). The proportional control gain is chosen so as not to saturate the actuator when the mass is at the extreme measurable limits of position. The static friction force bounds were based on the experimental estimates from Chapter 3. The sticking zone bounds are determined from these static friction levels and the proportional gain (3.4). The derivative and integral gain were chosen somewhat arbitrarily. The mass is assumed to be initially at rest outside of the sticking zone. 36 Table 4.1. Common parameters for all three PID simulations. Mass m 1 0 kg. PID K, 20.0 N/m Control K: 0.5 Ns/m. Gains Kt 40.0 N/ms. Static Friction F.“ 1.744 N. Bounds F.‘ -1.439 N. Sticking x, -0.0872 m. Limits x” 0.0719 m. Initial x0.) 0.15 m. Conditions V0,) 0.0 m/s. A single pair of equations can be used to define each of the three classes of slipping force models. These slipping force functions are represented as F;(v) = F; — AF+[1- (““91 + b‘v (4.9a) Fd‘(v) = F;-AF'rr-e""'5’]+b-v, (4.9b) where If and F: are the limiting values of the slipping force at zero velocity, AF " and AF ' are the respective drops from this level to the kinetic force level, and v: and v; are velocity constants defining the characteristic velocity at which 63% of the drop occurs. Where applicable, numerical values for the drop to kinetic friction are based on a kinetic coefficient of friction of 0.5, which roughly corresponds to a steel mass sliding on a wooden surface (Miller, 1977). The first PID time response uses the Coulomb plus viscous friction slipping force model (2.12). As seen from the parameters in Table 4.2, the slipping force linearly approaches the static friction bounds as velocity approaches zero from the left or right respectively. The artificial zero parameter is naturally chosen to be the assumed resolution of velocity measurementdue to sensor quantization. 37 Table 4.2. PID Simulation parameters for the Coulomb model. E.“ 1.744 N. F; (v) , 95+ 0.0 N. VI NA b+ 0.5 Ns/m. F.’ -1439 N. Far-(V) M?" 0.0 N. v; NA b’ 0.5 Ns/m. ArtificialZero ot 0.01m/s The simulated time response for PID control with the Coulomb plus viscous slipping force model is shown in Figure 4.3. The time scale is extended to 120 seconds since the mass is stuck most of the time. While stuck, the integral force term builds up until it reaches the respective static friction bound. At this point, the mass begins to move and gets stuck on the other side of the origin, as predicted by the theoretical analysis. x,x,,,x,: (x10“m.) F: (N.) ’4IIIITUIIIIITlIlTrIITTVUIIIITIUIUrrUI 0 20 40 60 80 100 120 Time (sec.) Figure 4.3. PID Control simulated time response with the Coulomb model. 38 From the displacement response and shape of the force response, the mass position appears to converge to the origin, though it looks like it will take a long time. As the error becomes smaller, the force due to the integral of the error requires more time to build up to a static friction level. From Figure 4.3, the time between successive “jumps” of the mass increases with each jump. This result indicates that PID control can be used to effectively regulate the position of the mass with a Coulomb plus viscous friction model, and will be investigated with additional detail in Section 4.3. The second PID time response uses the stiction plus viscous slipping force model. Similar to the previous slipping force model, and as indicated in (2.13), this model linearly approaches kinetic friction levels as the velocity approaches zero from the left or right respectively. The kinetic friction levels are arbiu'arily chosen for the simulation to be half of the static friction levels. The artificial zero parameter is again chosen to be the assumed resolution of velocity measurement due to sensor quantization. Table 4.3. PID Simulation parameters for the stiction model. P? 0.872 N. Flo) ép 0.0 N. v: NA b” 0.5 Ns/m. 12' -0.719 N. F;(V) ‘ AF' 0.0 N. v. NA b‘ 0.5 Ns/m. ArtificialZero ot 0.01m/s The simulated time response for PID control with the stiction plus viscous friction model is shown in Figure 4.4. The time scale is one tenth that of the previous plot (Figure 39 4.3), since the mass is no longer stuck most of the time. The mass still gets stuck twice per period of the apparent oscillation. 4 1 j .5": 2.. o! 3. 31‘ ,' F ‘. u i 1’} ." fl": 4 Il’: I r x" .1“ t 1414 11.: -1111 1.4 ._ 1. I 1'. | : ;: | I ll . l j: | t' ' i 0- ' ' — '. ' r.- ' ° - 1 | 1 l I 'J , '2 | 1, : I .4 '. J *W 1 .. _ _%’f— _ 134‘ _ 1‘ :1 ,' J .' )1 .5 ' ,5 ‘4 .I' I .3 l 'i f x U "" x, '2 : E. '...'\F :-.= x,x,,,x,: (x10'1m.) ‘ 17.12,: (N-) '4 I I I I I I I I I I I I I I I I I I I T I I r I I I I T I I TfT I I 0 2 4 6 8 10 12 Time (sec.) Figure 4.4. PID Control simulated time response with the stiction model. The overall friction force and the applied force is plotted in Figure 4.4. When the mass is stuck, the friction force is identical to the applied force until the static friction bound is reached. At this point, the mass begins to move and the slipping force drops to the kinetic level. The applied force overshoots the static friction bound thus increasing the acceleration of the mass (2. 15). Unlike the previous case, the position of the mass appears to settle into a stable limit cycle. This result indicates that PID control can not be used to effectively regulate the position of the mass with a stiction plus viscous friction model, and will also be further investigated in Section 4.3. The final PID time response uses an exponential plus viscous slipping force model (2.14). This model exponentially approaches the kinetic friction levels from the static level as the velocity increases away from zero from the left or right respectively. As in the 4o stiction model, the kinetic friction forces are chosen to be half the level of the respective static friction forces. Table 4.4. PID Simulation parameters for the exponential model. 1'? 1.744 N. F:(v) . Ap+ 0.872 N. v: 0.1 m/s. b” 0.5 Ns/m. F.‘ -1439 N. Fitv) . 95- -0719 N. V; -0.1 m/s. b’ 0.5 Ns/m. ArtificialZero ot 0.01m/s The simulated time response for PID control with the exponential plus viscous friction model is shown in Figure 4.5. As in the previous case, the mass appears to oscillate in a stable limit cycle motion and also gets stuck twice per period. 4 1 2-3 " 0'3 1 1 .2_: x,x,,,x,: (x10'1m.) : RF}: (N.) -4 rIIIWIIIIIIlIIIIIIIIIIIIIIIIIIIIrII 0 2 4 6 8 10 12 Time (sec.) Figure 4.5. PID Control simulated time response with the exponential model. 41 The overall friction force is again plotted with the applied force on this plot. Though the applied force overshoots the static friction bound when the mass breaks free, the friction force does not drop as rapidly as in the previous case, thus the resulting acceleration of the mass is not as large. These results also indicate that PID control can not be used to regulate the position of the mass with friction modelled by exponential plus viscous friction. 4.3 Limit Cycle Simulations Existence and stability properties of a limit cycle for an explicit slipping force function can be determined through the construction and examination of mappings of trajectories from points in S, to S, and then back to S, (Figure 4.2). In fact, existence and stability of limit cycle trajectories can be determined from mappings of the x-coordinate of a trajectory in S, into the subsequent x-coordinate in 5,. Simulation data of x-coordinates in S, mapped back into S, can be used to numerically construct a continuous iteration function for a given set of system parameters. Limit cycles will appear as fixed points of these iteration map functions, where the fixed point is the x amplitude of the cycle. The contraction mapping theorem can then be applied to prove the existence of fixed points (Vidyasagar, 1978, pp. 73-78). The equations for the stick-slip friction force (2.1-2.4), along with the three classes of models for E,P(v) (2.12-2.l4) completely define the friction force. Motivated by the time response simulations above it is useful to classify symmetric slipping force models with the following property (4.10). U = {Fw(v)| 3v. > 0 3 F;,,P(v) < E‘s/v e (0,v°)} (4.10) If a slipping force model belongs to the set U, it will take on values which are less than the static friction level in a region near zero velocity. Note that since the slipping force is assumed symmetric, we only need to look at positive velocities. It will be demonstrated 42 that models belonging to the set U will promote limit cycles in the 1-DOF mass system under PID control. Using the models defined in Chapter 2, the following classifications are obtained (Table 4.5). Table 4.5. Classifications for each of the three stick-slip models. Coulomb Model plus Viscous Damping pr(v) e U Stiction Model plus Viscous Damping Fd,’(v) e U Exponential Model plus Viscous damping pr(v) e U The slipping friction forces from the Coulomb model are always greater than or equal to the static friction levels therefore it is not in the set U. The stiction model has a discontinuous drop from the static level to a non-decreasing slipping force as the velocity increases and the exponential model undergoes a continuous drop from the static level to a lower kinetic level then increases, therefore both belong to the set U. Using the any of the three slipping force models, the system equations (4.2) are nonlinear. They are only piecewise linear for the Coulomb and stiction models. Since no general closed form solution exists for all three classes, a numerical simulation was performed to determine the contraction mappings. Three representative examples from each of the slipping force model classes were numerically simulated. All cases used the following set of common parameters (Table 4.6). 43 Table 4.6. Common parameters for PID limit cycle simulations. Mass m 1.0 kg. PID K, 10.0 N/m. Control K1 10.0 Ns/tn. Gains Kt 100.0 N/ms. Static and Kinetic F. 4.0 N. Friction Bounds FL 2.0 N. Artificial Zero 01 0.001 rn/s. Simulations were performed in double-precision on a Digital Equipment Corporation Micro-VAX, using a Hamming predictor-corrector general integration method with an initial step size of 0.001 s. For each run, a number of initial conditions inside S, were chosen to provide an evenly distributed data set for generating the iteration functions. Each x, is an x-coordinate in S, and x,“ is the x-coordinate for the next time the trajectory enters S,. A second order polynomial of the form: x“, = (a0 + a,x, + a,x,2), was fit through each set, and then plotted. There is one iteration function fit associated with each set of parameters. The numerical fixed points were obtained by guessing an initial value and then iterating the fitted function until the error between successive iterations was less than 1.0E-9. Convergence is guaranteed since the slopes of the iteration function over the fitted domain are less than unity. A conservative error estimate, 5, for each fixed point is obtained by using the largest absolute residual error from the curve fit to form a band around the fitted function. All data from the respective curve fit is therefore guaranteed to be inside this band. The fixed points of the two functions defining the edges of this band are then used to determine conservative error estimates in the original fixed point. These error estimates are presented with each of the fixed point estimates. For the Coulomb plus viscous friction models, the three representative cases were obtained by varying the viscous damping coefficient, b (T able 4.7). Table 4.7. Numerical results of fitted iteration functions for the Coulomb model. b a,(x10*) a, a,(x10") (x :1: 6)(x10'3) 2.0 1.87 0.814 6.37 1.01 :1: 1.51 1.0 2.05 0.853 7.39 1.40 i 2.10 0 0 1.99 0.906 7.00 2.12 :1: 4.65 0.5 ‘ 0.4- xr+r ‘ 0.3“ .1 0.2.4 .1 0.1— H b=0.0 HI b=1.0 H b=2.0 0. I I I I I I 0.0 0. 1 0.2 x, 0.3 0.4 0.5 Figure 4.6. Iteration functions for PID control with the Coulomb model. The simulation results for this class are plotted in Figure 4.6. Shown on this plot is the line: x”, = x, , along with the fitted iteration functions for the three test cases. A fixed point occurs wherever an iteration function crosses the x“, = x, line. All three test cases (Figure 4.6) have a fixed point at the origin, i.e. x“, = x, a: x, = 0. The magnitudes of the intercepts (a,) of the polynomials from the curve fit, as seen from Table 4.7, indicate that no limit cycles exist for this class of models. Within the accuracy of the simulation, all trajectories eventually converge to the origin (Table 4.7), and no limit cycle exists for the Coulomb friction model. 45 For the stiction plus viscous friction models, the three representative test cases were similarly obtained by varying the viscous damping coefficient, b (Table 4.8). Table 4.8. Numerical results of fitted iteration functions for the stiction model. b an a, 02 x :1: 5 2.0 0.0452 0.635 0.236 0.136 :1: 0.003 1.0 0.0489 0.663 0.246 0.165 :t: 0.004 0.0 0.0565 0.669 0.334 0.219 :1: 0.006 0.5 0.4-j xi+l “ 0.3— 0.2- 0.1— / / 13-6 b=l.0 - / x—x b=2.0 0.0 I I r l r I I 1 l 0.0 0.1 0.2 x,- 0.3 0.4 0.5 Figure 4.7. Iteration functions for PID control with the stiction model. The simulation results for the stiction plus viscous friction models are plotted in Figure 4.7. Each of these three test cases have distinct non-zero fixed points which are presented in Table 4.8 for each of the cases. From Figure 4.7, and as verified with the polynomial curve fit, the local slopes of these iteration functions near the fixed points are always less than unity. The contraction mapping theorem cantbe used to prove that these fixed points are all stable (Vidyasagar, 1978, pp. 73-78), thus stable limit cycles exist for the stiction plus viscous friction model. These limit cycles are reduced in magnitude with increasing viscous friction. 46 For the exponential plus viscous friction models, the three representative cases were obtained by varying both the viscous damping coefficient b, and the velocity constant v0 (Table 4.9). Table 4.9. Numerical results of fitted iteration functions for the exponential model. b V. do 01 02 x :1: 5 0.0 0.05 0.0628 0.732 0.190 0.297 :1: 0.008 0.0 0.10 0.0638 0.767 0.138 0.343 :t 0.004 1.0 0.10 0.0553 0.735 0.120 0.233 i 0.002 0.0 0.1 0.2 x,- 0.3 0.4 0.5 Figure 4.8. Iteration functions for PID conu'ol with the exponential model. The exponential plus viscous friction model results (Figure 4.8) also indicate the existence and stability of limit cycles but at larger amplitudes than those for stiction (Table 4.9). The local slopes near the corresponding fixed points are less than unity, so these fixed points are also stable, thus stable limit cycles will man for the exponential friction model. The Coulomb iteration function results indicate that as viscous damping is decreased, the slope of the iteration function increases and the fixed point remains near 47 zero. Similarly, the stiction iteration function results indicate that as the viscous damping is decreased, the fixed point of the iteration function increases. This observation is important since fast convergence to the fixed point (stable limit cycle) is obtained when the slope of the iteration function near a fixed point is close to zero. Since the cases investigated here all have slopes near unity, the convergence is slow. Raising the iteration functions further away from the x-axis (i.e. increasing 0,) has the effect of increasing the amplitude of the limit cycle. The presence of viscous damping can therefore help stabilize by either decreasing convergence times as in the case of Coulomb models, or by reducing the limit cycle amplitude as in the case of stiction and exponential models. 4.4 Experimental Investigation Using the experimental l-DOF mass system described in Chapter 3, the theoretical existence of limit cycle behavior was validated. As observed above, slipping force functions which do not have a definite drop from the static level will not admit a limit cycle response. A suitable braking material was found which produced a limit cycle behavior in the experimental system. This same brake was used in all the experimental investigations studied in this dissertation which required a brake. The original PD control algorithm with integral control added (Appendix B), was implemented on the same DEC LSI-l 1/23 computer, and ran at 200 Hz. For comparison with the PD response, all control parameters remained the same for the PID implementation except for the addition of integral gain. The pertinent control parameters for the experimental PID investigation are presented in Table 4.10. 48 Table 4.10. Experimental parameters for PID control. Mass m 1.375 kg. pro X. 20.0 N/m. Control K; 0.5 Ns/m. Gains K. 40.0 N/ms. Static Friction F.+ 1.744 N. Bounds F.’ -1.439 N. Sticking x, -0.0872 m. limits x, 0.0719 m. Initial 10.) 0.3 m. Conditions V0.) 0.0 m/s. The experimental system started up with the initial conditions given in Table 4.10. Under PID control, the mass moved toward the origin, and ended in a stable limit cycle. Due to the excessive amount of time required to reach this steady limit cycle, the time response for mass position and conu'ol force were not sampled from the initial startup conditions listed above. Voltage data was acquired fiom two channels with the use of An Apple Macintosh IIx computer and National Instruments LabView hardware and software at 1920 samples per channel and 8 Hz. per channel. The data was then scaled and plotted (Figure 4.9). 49 x: (x10'2m.) F: (N.) -3 IIIII'TITIIIIIIIT'IITIIIIIIITIITIII 0 40 80 120 160 200 240 Time (sec.) Figure 4.9. Experimental system time response for PID control. Several interesting observations can be made about the experimental time response of Figure 4.9. When the mass breaks free, the control force appears to be at the same level of the static friction forces which were estimated from the PD control experiment (Chapter 3). This agrees with the simulated predictions, and indicates that the estimates for static friction bounds are reasonable. As the static friction bounds are at different levels in the experiment, so too are the slipping forces. This can be seen by the difference in the amount of time the mass sticks on each side of the origin. For positive displacement error, the mass sticks for a longer time than for negative displacement error. The simulations were performed with symmetric slipping force functions except for different offsets, therefore the amount of time spent on either side of the origin was roughly equivalent. Meticulous observation of the experimental PID time response will prove to yield qualitative if not quantitative information about the character of the friction present. In this example (Figure 4.9), the static friction bounds can accurately be estimated, and we know that the slipping force is not of the Coulomb type. The slipping force must drop from the 50 static level to a lower kinetic level, and from the small amplitude of the limit cycle (approximately 1 cm.), we know this drop must be very small. Although this qualitative information does not provide a direct method for obtaining explicit parameters, it is still useful. Trial and error parameter tuning was used to produce a simulation in an attempt to match the above experimental result. The following slipping force parameters were found to produce a reasonable match (Table 4.11). Table 4.11. Slipping force parameters fitted to experimental limit cycle response. 1‘? 1.69 N. F;(v) Ap+ 0.0 N. v: NA 5' 20.0 Ns/m. F.‘ —1.39 N. Fir-(V) AF’ 0.0 N. v; NA b' 3.0 Ns/m. Artificial Zero 01 0.01 m/s 51 3 1 1 2': 2 - 1"". 53 3": J". Ii x I' ‘1. " " -’ l. i 5 ‘. 1 C. E 3 .‘, i f. a. f: k :' .l : 0' I I I ‘e . " Eu 3 3‘ ‘ '15 z I. i. I :0 a: 0’ :3 3 I. IV. A I 4" j a. f 4.1 :--.-1 ~ 4- — - 1 "- .3 ': r 1. i ’~ 1" ' “r— 0. ‘E— . 1 g: j: '. J I 1‘ 3 E .°. ’1': 2" ". ; i, : ( ' '1 .. 't|\ ' 1 " e r I 1 I I 2 F -2— _2 :J x: (x 10 m.) g F: (N .) -3 I I I I I I j—I T I I l I I I I I I I I I I I I I I I I I I I I I I 0 ' 40 80 120 160 200 240 Time (sec.) Figure 4.10. Simulated limit cycle response fitted to experiment. Initial conditions were chosen to match those of the experimental data. The simulated time response does not match up exactly with the experimental results, but the same character of the response is certainly present. The simple observations made above regarding the character of the slipping force model can be utilized to guide an engineer as to which models should be used. It is certainly important to use models which will actually admit the responses observed in an experimental situation. Chapter 5 Robust Nonlinear Compensation Putting two solids together is rather like turning Switzerland upside down and standing it on Austria - the area of intimate contact will be small. F.P. Bowden, B.B.C. Broadcast 1950. A simple but robust control law is still needed for the regulation of the 1-DOF system. PID control can only be used for some types of slipping force models, and therefore is not robust to small variations in stick-slip behavior. With PD control, there exists an undesirable set of multiple equilibrium points. This set can be reduced to a single point set by allowing the proportional gain to approach infinity, but this is not possible in practice. A discontinuous nonlinear compensation force is developed to supplement the PD control law. The resulting robust control law is provenito stabilize the 1-DOF system for any slipping force model satisfying the assumptions made in Chapter 2. This result is verified through numerical and experimental time responses. 5.1 Derivation of the Control Law The existence of a set of multiple equilibrium points, Em, (Figure 3.1) is undesirable, but the origin can be made into a unique equilibrium point of the system with the inclusion of a nonlinear compensation force 120:) in the PD control law F = -K;. - Kdv - EU), (5.1) where K P > 0 and K d > 0. Solving the l-DOF system equations (2.15) with the new control law (5.1) for the equilibrium points 52 53 i=0 => v=0 (5.2) v = 0 => [—K} - Kdv— F;(x)- Euph(v)- Ema—1w)» = 0 => -K,;t-F;(x)—Fm=0. (5.3) In order to have a unique equilibrium at the origin, we must choose a compensation force such that (5.3) is satisfied for 1:0 only. Choose the nonlinear compensation force to have the following form I»; = foc. (5. 4) Using the previous sticking limit definitions (3.4) with this compensation force and the equilibrium constraint (5.3), an equivalent requirement for uniqueness of the equilibrium is obtained by choosing an xc(x) such that the inequality xL S (x + xc) S x”, (5.5) is satisfied for x=0 only. This can be achieved with the following nonlinear function xc(x) O x>f,, (EH—x) O0. (5.7) xL=xL-e The addition of the nonlinear compensation force, defined by (5.4), (5.6), and (5.7), to the PD control law insures that the origin is a unique equilibrium point. The compensation force (5.4) is only active when the mass is between the extended sticking 54 limits EL 5 x S i”. The combination of the compensation force with the proportional force, K P(x + xc), can be viewed as a nonlinear proportional feedback control force (Figure 5.1). To simplify the implementation of this control force, the following two parameters are defined F“: = F: + er (5.8a) F,’ = F,‘ - K ,8. (5.8b) K; + F.(x> = K, o} and R' = {(x,v) 6 ms < o} (Siljak, 1969). Inside these regions, the right-hand side of (2.15) satisfies a 10cal Lipschitz condition, therefore we are guaranteed existence, uniqueness, and continuous dependence on initial conditions for solutions there (Vidyasagar, 1978, pp. 78-88). Trajectories are well defined and continuous until they hit the surface of discontinuity S where state derivatives are discontinuous. The notion of a “solution” of (2.15) must be generalized for trajectories which have discontinuous derivatives (Frlippov, 1964). This l-DOF system (2.15) and (2.1-2.4) has discontinuities of the first kind, or simple discontinuities. This means that there exists finite limits for the state derivatives as s —> 0+ and s —> 0‘ which do not coincide (Barbashin, 1970). It is easy to check from the limiting values of the state derivatives that at every point of S, the vector gradients of state 56 trajectories are directed with the same sense through the surface. Trajectories are never tangent to S and thus will always pass from one open region R*(R‘) through S into R’(R*). Trajectories are also never directed toward points in S from both open regions, therefore there cannot exist a sliding mode as found in the theory of variable structure systems. Over the entire phase space, trajectories are absolutely continuous. The qualitative behavior of trajectories in the phase space can be determined from the differential equations (2.15) using the bounds on the slipping force (2.4) and the nonlinear compensation force (5.4). Lemma: Every non-trivial trajectory will exit the quadrant it occupies in a clockwise fashion. The time spent in any one quadrant is bounded above by a constant determined by the initial condition. Proof: Consider solutions of the differential equations in each of the four quadrants of the phase space. In the first quadrant, x>O, v>O, and from (2.15), (2.4), (5.1), (5.6), Figure 2.2, and Figure 5.1 mv+ Kdv = -Kp(x+xC)-Fj(v) s E‘-b,v, (5.10a) and thus m9 + (K, + bo)v 5 if; < 0. (5.10b) From (5.10b) we know that v(t) is bounded above by a strictly decreasing function which approaches a negative value v(t) s {[v. — (E’) / O, and again from (2.15), (2.4), (5.1), (5.6), Figure 2.2, and Figure 5.1 mv+K,v = -K,(x+xC)-I«;+(v) 2 i;*-F,,*—b,v, (5.12a) and thus m} + (K, + b,)v 2 (if,+ -- If) 2 (i3;+ - Ff) = er > o. (5.12b) Similarly, from (5.12b) we know that v(t) is bounded below by a monotonic function which approaches a positive value v(t) 2 {[v - (er) / (K, + b,)]e""‘”’""”‘ + (er) / (K, + b,)}. (5.13) Trajectories must also leave the second quadrant in a clockwise fashion. Proof of this lemma in quadrants III and IV proceeds as in quadrants I and II respectively. Q.E.D. 5.3 Stability Theorem and Proof Now that a control force has been determined such that the l-DOF system (2.15) has a single equilibrium point at the origin of the state space, we need to show that it is stable. Theorem: The origin of the state space for the system given by equations (2.15) and 2.1-2.4), with the control force given by equations (5.1), (5.4), (5.6), and (5.7), is a globally asymptotically stable equilibrium point. Proof: Choose the following Lyapunov function candidate, which closely resembles the actual energy function of the system 58 V(x,v) = yzprz + va2 + g(x), (5.14) where KK,(ig)2 x > i” Kp(i,,x-}5x2) 0 51: Si” Kp(i,_x—}éx2) iLSxSO' )4K,,(5c',)2 x < i, 8(X) = 4 (5-15) l gov) 72nd,)” I | | | l l 1 5? <— H.)—————— I“ Figure 5.2. Nonlinear addendum, g(x), to the Lyapunov function candidate. It is easy to check that this Lyapunov function candidate is positive definite, decrescent, and satisfies a global Lipschitz condition. This function is also continuous with respect to xandv. We must now evaluate the time derivative of the Lyapunov function candidate, V(x,v), along solution trajectories of (2.15). The time derivative of V is the dot product of the gradient of V with the vector of state derivatives. As mentioned above, the state derivatives are discontinuous at all points in S, and from (5.15) the gradient of V is discontinuous at x=0. The derivative V(x,v) is discontinuous or does not exist for all points of S, however it exists and is well-defined for all points in the open regions R+ and R'. 59 A modification of Lyapunov’s direct method must be used for non-differentiable Lyapunov function candidates, such as the one presented above (Solncev, 1951; Hahn, 1963). The global asymptotic stability theorem of Lyapunov for continuous systems is modified by replacing the derivative V(.) with the Dini-derivate D‘V(.), where the * represents any of the four possible Dini-derivates (Rouche, et al., 1977). At any point where the derivative V(.) exists, all four Dini-derivates will have a common value equal to the derivative at that point (McShane, 1947). D'V(.) = V(.) V(x,v) e s. (5.16) Along solution trajectories, excluding points in S, we can compute the time derivative of the Lyapunov function candidate as V(x,v) = ng + v(F - F,) + g’(x)v = K gv + v(-K;t - K ,v — K ,Jc - Euva) - Emu - l(v))) + g’(x)v = - ,v2 — thup(v) + v(g’(x) - K rte), (5.17) where ' 0 x > i” Kp(i,,-x) O 0 where, along some trajectory (x,v), 1imV(t) = V, (5.20) [-50- The level set V(x,v) = V, is a closed curve in the phase space. From the previous lemma, trajectories must spiral clockwise around this curve, and must traverse the cycle within a fixed period. Along solution trajectories, going from to to t1 , V|‘ .-_- ("th = —f'(K,v2+ vii») dt, (5.21) and from (2.4) we have, VFW 2 bov2 Vv, therefore v|:‘ s —(K, + b,)j"'v2 d1: = -(K, + b,)j:v dx. (5.22) The last integral in (5.22) is just the area contained inside the path from x0 to x, in the phase space. Each limit cycle will decrease V(t) by at least an amount which is proportional to the area contained within the limit cycle, therefore no such V, can exist. Q.E.D. 61 5.4 Simulation Results A numerical simulation of the nonlinear stick-slip friction compensation force applied to the mass system will verify the stability of the origin. As in the previous chapter, in order to obtain a more realistic simulation, quantized states (simulating the measurement process) are used to compute the control force, which is in turn quantized (simulating the output process). In order to numerically verify robustness of the control law, a representative slipping force model from each of the three classes will be tested with the same control law. All simulations were performed in double-precision on a Digital Equipment Corporation Micro-VAX, using a Hamming predictor-corrector general integration method with an initial step size of 0.001 s. (Appendix A). Using similar assumptions as those in Chapter 4, the following parameters were chosen to approximately model an experimental system (Table 5.1), and were common to each of the three numerical simulations. 62 Table 5.1. Common parameters for all nonlinear compensation simulations. Mass m 1.0 kg. PD Control K, 20.0 N/m. Gains K4 0.5 Ns/m. Static Friction F.” 4.2 N. Bounds I? -4.0 N. Sticking xL -0.21 m. Limits IL 0.20 m. Nonlinear 8 0.005 m. Compensation 1E: 4.3 N. Force 13:- -4.1 N. Extended 55,, -0.215 m. Sticking Limits 55L, 0.205 m. Initial 2:0.) 0.3 m. Conditions V(ta) 0.0 m/s. Theoretically, any positive value for 8 can be used for the nonlinear compensation force. In a real implementation, 8 must be chosen carefully. This control law assumes that the discontinuous jump in the input force can actually be produced instantaneously. In reality, the bandwidth of the actuator will limit the bang-bang nature of the control. In (5.8), we have seen that K p8 is a force which is added to the static friction levels to insure that the net force on the mass is non-zero in the sticking zone (Figure 5.1). In this example, the static friction force bounds are known exactly. In actuality, there will only be estimates for these forces. Values for 13;” and if should be chosen to exceed the upper bounds of the estimates for the static friction force levels. The value of e in Table 5.1 corresponds to assumed measurements of the static friction levels with less than a 2.5% error, thus the control force levels, if and E", are chosen to be 0.1 N above the static friction levels. Again, the mass is assumed to initially be at rest outside the sticking zone. The slipping friction force function used in this numerical study is the same general representation as used in Chapter 4 (4.9) since each of the three classes of slipping force 63 models are easily represented by these functions. The first model is the Coulomb plus viscous friction model. The defining parameters for this model are given in Table 5.2. Table 5.2. Simulation parameters for nonlinear compensation with the Coulomb model. F.” 4.2 N. FR") ‘ éE* 0.0 N. v: NA b+ 0.5 Ns/m. F.” .4.0 N. Fifi) AP" 0.0 N. v; NA b' 0.5 Ns/m. ArtificialZero a 0.01m/s A simulated time response verifies the stability of the nonlinear compensator for the Coulomb plus viscous slipping force model (Figure 5.3). Several important observations can be made regarding this time response. As soon as the mass enters the extended sticking zone, the control force takes on its bang-bang character. Note that since the velocity is small, the derivative term in the control force is small also. M U) H l H I 03 I I ,--.-.--..-------..-:J x1. 0 UI llLlllJlLlllllllllll lllllllJIJJUI I ,.-"\ XE x,i,,,J'EL: (x10'1m.) F F, 15;: (N.) TTIIIIIIFIIIIIIIIIIIIIIrIIIIIII 0.5 1.0 1.5 2.0 Time (sec.) #1 P c Figure 5.3. Simulated time response for nonlinear compensation with Coulomb model. Whenever the velocity becomes zero, there exists a possibility of becoming stuck. From Figure 5.3, we see that the control force always exceeds the friction force when the velocity becomes zero, thus the mass continues moving. Finally, notice that the mass displacement does not overshoot the desired reference by very much. This is a characteristic of systems with a Coulomb type of slipping force. The second model is the stiction plus viscous fiiction force. The defining parameters for this model are given in Table 5.3. 65 Table 5.3. Simulation parameters for nonlinear compensation with the stiction model. 1'? 2.1 N. Fd+(v) ‘ AF" 0.0 N. VI NA b“ 0.5 Ns/m. 1'? -2.o N. Fd-(V) AF‘ 0.0 N. v: NA b' 0.5 Ns/m. ArtificialZeto a 0.01 m/s The time response for the 1-DOF system with a stiction plus viscous slipping force model and the nonlinear compensator is shown in Figure 5.4. Again, when the mass enters the extended sticking zone, the control force takes on its bang-bang character. For this slipping force model, both the displacement and force oscillate more than in the previous example before stabilizing, thus more control energy is required and it takes slightly longer to reach the desired equilibrium. ? 1?:::tf"‘ L-..— -- 5 ___.L ‘— I \. T- XL x,i”,iL: (xlO’lm) 17.151 (NJ I I I I I I I l I I I I I I I I I I I I I I I T I I I r I I 0.0 0.5 1.0 1.5 Time (sec.) 2.0 Figure 5.4. Simulated time response for nonlinear compensation with stiction model. Again, the control force always exceeds the friction force when the velocity becomes zero, thus the mass continues moving. For this model, the mass displacement does overshoot the desired reference and oscillates about the reference until stabilizing. The final model is an exponential plus viscous friction force, and is more similar in form to the stiction model than the Coulomb model. The defining parameters for this model are given in Table 5.4. Table 5.4. Simulation parameters for nonlinear compensation with the exponential model. 1‘? ‘ 4.2 N. F;(v) , Ap+ 2.1 N. VI 0.1 m/s b+ 0.5 Ns/m. FL‘ 4.0 N. F;(v) , 915- .2.0 N. V; -0.1 m/s b’ 0.5 Ns/m. Aru'ficialZero a 00le The time response for the l-DOF system with an exponential plus viscous slipping force model and the nonlinear compensator is shown in Figure 5.5. This response, as expected, is more similar in character to the stiction model response than the Coulomb model response. ’5‘: ”gr\ JaimiL: (x10'1m.) : F 17.15;: (N-) -7 I T I I I I I I T I I I I I I l T I j I I I l I I I I I I T 0.0 0.5 1.0 1.5 2.0 Time (sec.) Figure 5.5. Simulated time response for nonlinear compensation with exponential model. Again, when the mass enters the extended sticking zone, the control force takes on its bang-bang character. The displacement and force oscillate more than in the Coulomb example, but not as much as in the stiction example. Less control energy is required for stabilization, and it takes about the same amount of time to reach the desired equilibrium as in the case of stiction. Each of the simulation plots have a common characteristic. The bang-ban g effect of the control force is one where the force switches sign whenever the mass passes through the desired reference. For example, the force may be positive for a period of time until the mass overshoots the origin, then it becomes negative for a shorter period of time. This decrease in the time width of the pulse is observed in each successive pulse for each of the three examples. From the simulation results, the desired reference position of the 1-DOF system is a stable and unique equilibrium point with the nonlinear compensation force supplementing PD control. As long as reasonable estimates for the upper bounds on the maximum static 68 friction forces can be obtained,'this control technique is robust to the exact character of the slipping force. No limit cycle behavior is observed as in the case of PID control with the stiction and exponential models. 5.5 Experimental Results The nonlinear proportional feedback control law (5.4), (5.6), and (5.7), has also been tested on the experimental system to verify the stability of the origin. The estimates for the static friction bounds from Chapter 3 (Table 3.1) are used to determine the compensation force levels. In an experimental implementation of this compensation force, since we only have estimates for the static friction bounds, we do not need to choose an e explicitly. in order to insure that no sticking would occur in the sticking zone, the applied nonlinear control forces, I? and F" , only need to be chosen to exceed the largest estimate. The proportional and derivative gains were chosen identical to the gains in the simulation section for direct comparison. The initial conditions were also chosen identical to those from the simulations. 69 Table 5.5. Experimental parameters for nonlinear compensation. Mass m 1.375 kg. PD Control K, 20.0 N/m. Gains K; 0.5 Ns/m. Static Friction P? 1.744 N. Bounds P;- -l.439 N. Sticking XL -0.0872 m. Limits x41 0.0719 m. Nonlinear .+ 2.0 N. Compensation i: -2.0 N. Extended 3,, -0.1 m. Sticking Limits &, 0.1 m. Initial 10.) 0.3 m. Conditions V0.) 0.0 m/s. The velocity state, which is needed to compute the derivative control force, was estimated with a finite difference algorithm using position data. The floating—point controller algorithm was implemented on a DEC LSI-11/23 computer in Fortran IV, and run at 200 Hz. In the implementation, a first order digital filter was used in the estimation of the velocity state, but no filtering was used on the position data. The experimental data was acquired with the use of an Apple Macintosh IIx computer and National Instruments LabView hardware and software. From each channel, 2000 samples were taken at 1000 Hz. per channel. To experimentally verify stability and robustness of the nonlinear compensation technique, two tests were performed. In the first test, the friction brake was applied, since the estimates for static friction were obtained under this condition. It was observed in Chapter 3 that PD control led to a non-zero steady-state error, and in Chapter 4, PID control led to a limit cycle response when the friction brake was applied. The resulting time 70 response for this experiment is shown in Figure 5.6. Note the decrease in the width of the control force pulses as observed in the simulations. 4. I x .. 2 x11 2: 1” § 0: — — — -— —- ‘..: w- I xL -4; 3 F -6: ~ ~. -1 _ x,x,,,x,_. (x10 m.) I F: (N.) -8 I I I I I I l l' Ij I I I I I l I I I T I I l I I I I I I I 0.0 0.5 1.0 1.5 2.0 Time (sec.) Figure 5.6. Experimental time response for nonlinear compensation with braking. The time response of the experimental system with braking (Figure 5.6) is similar in character to the simulation results for a Coulomb plus viscous friction model (Figure 5.3). This result agrees well with the results in Chapter 4 obtained through the PID simulations and experiment. It was observed that the friction present in the braked system was stiction with a very small but definite drop from the static levels. The drop was estimated to be so small as to be almost Coulomb friction, but the definite drop was required to give a limit cycle response for PID control. The friction brake was completely removed for the second experiment. Because of this, the static friction bounds obviously decreased, but the remaining stick-slip friction was observed to be much less uniform with position. Theoretically the controller is robust enough to handle thisnew situation with no change. The resulting time response for this experiment is shown in Figure 5.7. 71 4 -1 E 55,, 21 \ OZ 1--------------------- ------ ----1 -2-5 1‘ 2 EL -43 6.? F I x3052”: (x10'1m) 3 F: (N.) -8 I j I I I rTTT I I I I I I l I T I I I I I l I I I I I I I 0.0 0.5 1.0 1.5 2.0 Time (sec.) Figure 5.7. Experimental time response for nonlinear compensation without braking. The resulting unbraked time response verifies the stability of the origin. Due to the overall decrease in the friction present, the mass reaches the origin much faster than in the braked situation. Chapter 6 Conclusions 6.1 Summary The physical phenomena of the stick-slip motion of one surface over another has been investigated. The friction force induced by the surface contact has been identified as the causal agent, thus it is termed a stick-slip friction force. A good mathematical model for this stick-slip friction force is presented which, when incorporated into a dynamic system model, will provide the stick-slip dynamics observed in practice. This model is composed of a sticking force and a slipping force which do not act simultaneously. The sticking force model provides values of the friction force when the relative velocity of the surfaces is zero, and the slipping force model is assumed to be a function of the relative velocity. The common models for friction exclusively neglect the sticking force altogether for reasons of analytical simplicity. Omission of the sticking force from the overall friction model will certainly reduce the richness of the set of possible dynamics obtainable from the model. Three of the common slipping force models and their characteristics are investigated throughout this analysis to determine and compare their effects with different control strategies on the regulation of a single degree-of-freedom system. These three models were Coulomb, stiction, and exponential. All three included a linear viscous damping term since it is the most common friction model. Using a Proportional+Derivative (PD) control law, there exists a set of multiple equilibrium points around the desired equilibrium. This result is verified theoretically as well as experimentally. The steady-state error bounds are inversely related to the 72 73 proportional gain. Since this gain must be finite, there will always exist non-zero steady- state error bounds for PD control. The static friction bounds for the stick-slip friction force in a given l-DOF system can be experimentally estimated using the PD control law. With each new proportional gain, there are two new steady-state error bounds which can easily be measured. From these measured bounds, the static friction bounds can be computed As long as the static friction forces are not dependent on position, the estimates obtained through this procedure will be very good. Each of the estimated static friction bounds for the experimental system with the friction brake applied were no more than 3% away from the mean values of the respective samples. This indicates that there is only a minimal dependence on position, and the estimates are very reasonable. A Proportional+Integral+Derivative (PID) control law was investigated in an attempt to remove the steady-state error. The initial theoretical analysis demonstrated the existence of trajectories which appeared to have a limit cycle behavior. The nonlinearity of the stick-slip friction force did not allow a closed form solution, therefore numerical simulations were utilized to provide the limit cycle existence results. Simulations with the stiction and exponential models appeared to exhibit limit cycle behavior over the finite time interval of the simulations. The Coulomb model did not appear to exhibit a limit cycle behavior. These results were verified through the numerical construction of a mapping of points in the sticking zone, where the fixed point of this mapping is the limit cycle amplitude. Several representative cases from each class of slipping force model were tested under the same conditions. A Coulomb friction model in the l-DOF system with PD) control does not exhibit a limit cycle response, whereas the stiction and exponential models promote definite limit cycle dynamics. This fact can be used advantageously to determine the qualitative character of the slipping force present in a real experimental l-DOF system. If a limit cycle exists 74 under PID control, the slipping force cannot have the characteristics of the Coulomb model. The amplitude of the observed limit cycle depends somewhat on the amount of the drop from the static friction level to the kinetic friction level. As the drop decreases to zero, the limit cycle amplitude also decreases to zero. A discontinuous nonlinear compensation force is developed to supplement the PD control law. This compensation force is constructed to make the origin a unique equilibrium point of the system. The discontinuous compensation, combined with the natural discontinuity of the stick-slip friction force create a surface of discontinuity in the phase space. A sliding mode does not exist, and since the discontinuities are simple, all trajectories are absolutely continuous. The origin of the 1-DOF state space is proven to be globally asymptotically stable by two Lyapunov methods. The first uses a modification of the direct method for discontinuous systems. Dini-derivates must be used instead of derivatives for non-differentiable Lyapunov function candidates. The second method shows that the energy of the system (given by the Lyapunov function) must decrease a finite amount with each cycle about the origin, therefore no sustained limit cycle can exist. Stability of the origin is verified through numerical simulations and an experimental demonstration. Each of the three classes of slipping force models were used in a simulation study with the same nonlinear control law, and the same operating conditions. All three test cases satisfied the assumptions that the slipping force functions be Lipschitz in their respective domains, and they lie within prescribed bounds. Stability was verified with each of the three models, thus the nonlinear compensation is robust to the character of the slipping force. Similar to the results with PID control response, the qualitative character of an experimental slipping force can be determined by comparison with the simulated responses for each of the three classes. The braked and un-braked experimental responses verify the stability of the nonlinear compensation force as well as its robustness. 75 6.2 Recommendations for future work The robust nonlinear compensation technique works well for the 1-DOF system, but its effectiveness on a multi-DOF system is questionable. Preliminary investigations of the extension of this control law applied to a 2-DOF system (the inverted pendulum) with stick- slip friction present indicate that it cannot be simply extended. This can be explained by looking at the sticking force. When one of the degrees of freedom is “stuck” the nonlinear control law injects enough energy into the system to get it un-stuck. This can adversely affect the resulting motion as indicated by the simulation study with the inverted pendulum. It is probable that a greater degree of success is attainable with this compensation technique for systems where an actuator and sensor pair are directly associated with each degree of freedom. These types of systems are commonly found in the robotics industry. A natural extension of this work is to develop online or offline observation algorithms for identification of the slipping forces and static friction levels in an experimental system. The qualitative characterizations which can be obtained through the PID time responses or the nonlinear compensation time responses will provide a better “first guess” to the exact character of the slipping force present. The nonlinearity associated with real stick-slip friction tends to be very complex and . dependent on many quantities. The static friction forces can be strongly dependent on position, and even the slipping forces may vary somewhat with position. With this in mind, successful observation and identification may only be attainable with the “nicest” experimental systems where the slipping forces and static friction bounds are uniform. In dealing with real systems, it is better to investigate the control and observation techniques which do not attempt to identify the exact friction present, but are robust enough to successfully stabilize the system. Appendices Appendix A Simulation Models A.1 Numerical Simulation Numerical solutions to the non-linear differential system equations were computed on a MicroVAX-II running VMS V5.2 using the simulation package DIFFEQ (CCUG, 1989). All computations were performed in Double Precision. The following FORTRAN- 7 7 Functions/Subroutines were linked into the DIFFEQ program for the actual simulation. A.2 Stick-Slip Friction Model. The following function computes values of the sticking force when the velocity is within the artificial zero range as defined in Chapter 2, Section 2.2. When the velocity is outside of the zero range, the slipping force is evaluated. C C .................................................................. c ----------------------------- Fd (V) ------------------------------ C FUNCTION Fd (V, Eps, Fsp, an, Fr, Cp, Cn) C C This function computes the value of the damping force. C C ---Variable Identification C C V the sliding velocity between surfaces C Eps artificial zero ( > O ) for computation purposes C Fsp Positive ( > 0 ) static friction force C an Negative ( < 0 ) static friction force C Fr Reaction force (see notes for definition) C Cp coefficient array for pos. non-zero function C Cn coefficient array for neg. non-zero function C C ---Declare variables C REAL*8 Fd,V,Eps,Fsp,an,Fr,Cp(4),Cn(4) C C ---Compute the damping force C IF (V .GT. Eps) THEN C C ---Velocity is positive 76 77 C Fd = FlV-Eps,Cp) C ELSE IF (V .LT. -Eps) THEN C C ---Velocity is negative C Fd = F(V+Eps.Cn) C ELSE C C ---Velocity is zero C IF (Fr .62. Fsp) THEN C Fd = Fsp C ELSE IF (Fr .LE. an) THEN C Fd = an C ELSE C Fd = Fr C ENDIF C ENDIF C RETURN END C C _____________________________________________________ - c ----------------------------- Fd(V) --- C The slipping force is computed with the following function whenever the velocity is non-zero, or outside of the artificial zero range. The value of the slipping force is a function of the velocity, and a set of four parameters stored in a parameter vector. With a proper choice of the parameters, any of the three general stick-slip friction models presented in Chapter 2, Section 2.3, can be represented. This function can produce Coulomb, Stiction, or Exponential models with viscous friction added in. C C ................... -. —————————————————————————————————— c ----------- - ------ F(V,C) ------------- C FUNCTION F (V, C) C C This function computes the value of the damping force for C nonzero velocity. C C --—Declare variables C REAL*8 F.V,C(4) O 78 C ---Compute the damping force C F = C(1) + C(2) * (1.000 - DEXP( -V / C(4))) + C(3) * V C RETURN END C C ———————————————————— ——— — - c ----------------------------- F(V,C) -- --------- C A.3 One-Degree-Of-Freedom Model The l-DOF system equations were written into the following subroutine for evaluation by the DIFFEQ package. Through this package, a set of parameters could be interactively changed to obtain different modes of operation. For example, the PD and PID controllers were both included in the following code. 0 0 0 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 'U N U) H I! I I I I I I I I I I I I I I I I I SUBROUTINE PZSIM This subroutine has the simple forced mass problem with the stick-slip damping model, controlled by PID feedback. ---DIFFEQ subroutine name 00000000 SUBROUTINE FCT ( T, 5, Bars ) C C -—-Variable Identification C C T The real time variable C S State vector C DerS Derivative of state vector C PA Parameter vector C C Kp Position control gain C Kd Velocity control gain C Ki Integral control gain C M Actual mass C C Fsp Positive static friction force C an Negative static friction force C Cp Coefficient array for damping force (+) C Cn Coefficient array for damping force (-) C dF Damping force from the model C Fd Name of function generating damping force C C Fr Reaction force (or control force) C Fc Compensation force C Fpid PID control force C C Est Quantization precision for position C EpsV Quantization precision for velocity C EpsF Quantization precision for force 000000000000 0 0 0 0 0 0 0 EpsT Tk X V Y E DE IE 79 Quantization precision for time Quantized time variable Quantized mass position Quantized mass velocity Quantized integral of the error Position error Velocity error Integral of the error ---Declare variables REAL‘B REAL‘B + REAL'B REAL*8 ---Common COMMON DerSIZO),S(20),PA(20),T Kp,Kd,Ki,M,Fsp,an,Cp(4),Cn(4), Fd,dF,Fr,X,V,Y,E,DE,IE,Xd,Vd, Est,EpsV,EpsF,EpsT Fpid,Fc,Xc,Xh,Xl,Eps,G,dG Tk,Told,Xold,Vold for fct parameters /FCTCOM/PA ---Set up all local parameters EQUIVALENCE (Kp, PAIIII: + + + + + + + + + + + (Kd, PA(2)). (Ki, PA(3)). (Est, PA(4)). (EpsV, PA(5)). (EpsF. PA(6)). (EpsT, PA(7)). (Eps, PA(8)). (Fsp, PA(11)). (an. PA(12)I. (Cle),PA(13)). (Cn(l),PA(17)) —--Initialize local constants DATA M,Xd,Vd /1.0D0,0.0D0,0.0DO/ ---Simulate measurement by quantizing the states Tk a Q ( EpsT, T ) X = Q ( Epsx, 5(1) ) V Y Q ( EpsV, 5(2) ) 5(3) --—Compute the errors E = Xd - X DE = Vd - V IE = - Y ---Compute the PID control force Fpid - Kp*E + Kd*DE + Ki*IE ---Compute the compensation force 80 Xh = - (an / Kp) + Eps X1 = - (Fsp / Kp) - Eps C IF (X .LT. X1) THEN Xc - 0.000 ELSE IF (X .GT. Xh) THEN Xc = 0.000 ELSE IF (X .GT. Est/2.0DO) THEN Xc = Xh - X ELSE IF (X .LT. -Est/2.000) THEN Xc - Xl - X ELSE Xc = 0.000 ENDIF C PC = - Kp * Xc C Fr = Q ( EpsF, (Fpid + PC) ) C C ---Compute the actual damping force C dF = Fd (5(2), EpsV/2.0DO, Fsp, an, Fr, Cp, Cn) C C ---Eva1uate system derivatives C IF ( ABS(5(2)) .LE. EpsV/2.000 ) THEN 0er5(l) = 0.000 ELSE 0erS(1) = 5(2) ENDIF C Der5(2) = (Fr - dF) / M C 0erS(3) = - E C C ---Compute auxiliary outputs C 5(4) = 10.0 * 5(1) 5(5) = 10.0 * 5(2) 5(6) = Y 5(7) = Fr 5(8) = dF 5(10) = Fc C C ---Compute the Lyapunov function C IF (5(1) .GT. Xh) THEN G = Kp * Xh * Xh / 2.000 d6 = 0.000 ELSE IF (5(1) .LT. Xl) THEN G = Kp * X1 * X1 / 2.000 d6 = 0.000 ELSE IF (5(1) .GE. 0.000) THEN G = Kp * 5(1) * (Xh - 5(1)/2.000) dG - Kp * (Xh - 5(1)) ELSE G = Kp * 5(1) * (X1 - 5(1)/2.000) d6 = Kp * (Xl - 5(1)) ENDIF C 3(11) = G + ( Kp*S(l)*S(l) + M*S(2)*S(2) I / 2.0D0 5(12) = Kp*$(1)*S(2) + M*S(2)*Der$(2) + dG*S(2) 5(13) = 10.0 * Xh 5(14) = 10.0 * Xl 81 C RETURN END C C .................................................................. C __________________________________________________________________ C C ---Insert all the INCLUDEd functions C INCLUDE '[SOUTHWARO.RE5RCH.PID.COMMON]FD.FOR' INCLUDE '[SOUTHWARO.RESRCH.PID.COMMON]F.FOR' INCLUDE '[SOUTHWARD.RE5RCH.PID.COMMON]Q.FOR' C C ——————————— - - ....... C .......................................... C In the FCT subroutine above, the actual position and velocity are quantized to simulate the discrete measurement process. This quantization is performed by the following function . C C ................................. - - .. - ........... C ---------------------------- Q(Eps,X) ---------------------------- C FUNCTION Q ( Eps, X ) C C This quantizer maps real X values into the Q—Eps space. C C ---Declare Variables C REAL*8 Q,X,Eps,MAX C C ---Compute the maximum value quantizable C MAX = Eps ‘ (2.000 ** 30) C C ---Compute the R value C IF (ABS(X) .LT. MAX) THEN Q = Eps * JIDINT( (X/Eps) + 05IGN(0.500,X) ) ELSE Q = X ENDIF C RETURN END C C ____________________ — - - — ...... C ---------------------------- Q(Eps,X) ---------------------------- Appendix B Experimental Controller 8.] Controller Implementation The two control algorithms used in this investigation, PID and PD with nonlinear compensation, were implemented on a DEC PDP-11/23 running RSX-11M+. The programs were written in FORTRAN IV and MACRO-11. All control computations were carried out in Single Precision Floating Point. Four modules were common to both algorithm implementations. These four modules are IPOSN, PAUSE, SETICK, and SAMPLE, and can be found in Section B3. The specific PID routines are in the following section. 8.2 PID Controller. The PID control algorithm is implemented in the following program. A single input voltage is read and filtered with a first order digital filter. The derivative of this filtered signal is computed and filtered as well at user specified time constants. The integral of the input error is also computed. The PID control signal is then computed and applied. This program allows the user to interactively tune all parameters. c C —————————————————————————————————————————— c --------------------------------- PID - --- ----- c PROGRAM PID This is the main program which accesses all the control subroutines for PID control. Written by: STEVE C. SOUTHWARD ---Program Calling Tree PID IPOSN 0000000000 82 000000000 0 0 0 00 SAMPLE CGAINS CONTRL RUNLSF PAUSE SETICK ---Declare Variables REAL G(10).F(10) LOGICAL*1 IDATA DIMENSION IDLIST(21) EXTERNAL SETICK 83 COMMON /ADDR1/ IKSR,IKDR,INCH,ITCK COMMON /ADDR2/ IPSR,IPDR,IBPR,IADVEC,IADERR COMMON /ADDR3/ ICSR,IADR,IDAR,IADBUF COMMON /MISC2/ F,G,IZREF COMMON /TIME/ I.IEND ---Initialize the zero reference DATA IZREF/O/ ---Initial PID control gains, and other parameters F(l) - 10.0 F(Z) = 5.0 F(3) - 2.0 H F(4) ‘ .0/(1.0 + 3.0) F(S) = 1.0 - F(Q) F(6) = 1.0/(1.0 + 10.0) F(7) = .0 - F(6) I—I F(8) = -5.982 F(9) 1.68 F(10) 8 1.0 ---Set up address registers IPSR - "102 IKSR = "177560 IKDR = IKSR + 2 IPSR - IKSR + 4 IPDR ‘ IKSR + 6 ICSR = "170420 IBPR ‘ ICSR + 2 ICKV s "440 IDAR = "170440 IADR = "170400 IADBUF - IADR + 2 IADVEC - “400 IADERR = IADVEC + 4 IBELL - "7 IZERO = "4000 ---Protect System from crash !Kp \ !Kd PID controller gains !Ki / !position !filter !derivative !filter !input scale factor !output scale factor !maximum allowable error !processor status word !keyboard status register !keyboard data register !printer status register !printer data register !real time clock status reg. !clock buffer preset reg. !clock vector address !D/A base output addr.: ch.0 !A/D control status register !A/D buffer register !A/D done interrupt vector !A/D error interrupt vector !bell tone !zero value for D/A converter 0 0 84 IDLIST(1) - IKSR !keyboard status IDLIST(2) - IPEEK(IKSR) IDLIST(B) - IPSR !printer status IDLIST(4) 8 IPEEK(IP5R) IDLIST(S) = IDAR !zero D/A channel 0 IDLIST(G) - IZERO IDLIST(7) - IDAR + 2 !zero 0/A channel 1 IDLIST(B) ' IZERO IDLIST(9) - IDAR + 4 !zero D/A channel 2 IDLIST(IOI ' IZERO IDLIST(11) = IDAR + 6 !zero D/A channel 3 IDLIST(12) - IZERO IDLIST(13) - IADR !A/D status IDLIST(ld) - 0 IDLIST(IS) - IADVEC !A/D done interrupt vector IDLIST(16) = IPEEK(IADVEC) IDLIST(17) - ICSR !clock status IDLIST(18) ‘ IPEEK(ICSR) IDLIST(19) 8 IPSW !processor status IDLIST(20) - IPEEK(IP5W) IDLIST(21) = 0 !end of list marker CALL DEVICE(IDLIST) ---Format Statements 500 FORMAT(/,/,/,2X,'Discrete Time PID Control', + /,2x,' ------------------------- ') 510 FORMAT(/,2X,'Enter the INPUT A/D channel: ',5) 530 FORMAT(/,2X,'Enter the OUTPUT D/A channel: ',$) 540 FORMAT(IZ) 550 FORMAT(/./,10X,' ----- MAIN MENU ----- ', /,/,2X,'(I)nitialize the zero reference position', /,/,2X,'(G)ain settings for PID control I filters', /./.2X,'(5)tatus of Zero position', /,/,2X,'(R)un the active PID controller', /,/,2X,'(Q)uit the program', /,/,2X,'Choose one of the above. . . ',S) 560 FORMAT(/,/,2X,'5TATUS of Zero Reference. . .', + /,2x,' ------------------------------ ') 570 FORMAT(/,2X,'INPUT zero reference: ',I6) 590 FORMAT(/,/,2X,'Exiting Control Routine. . .',5(/)) 600 FORMAT(/,/,2X,'Zero reference position NOT initialized') + + + + + + ---Calibrate the PAUSE subroutine IEND = 10000 !set initial value IPR = 7 !highest priority interrupt 10 I = INTSET(ICKV,IPR,1,5ETICK) !attach to RTC vector IF (I .EQ. 0) GOTO 20 WRITE(7,*)'INT5ET error -- RTC vector, CODE - ',I GOTO 10 20 IRATE = "111 !set up clock to generate ICOUNT - -16667 !an interrupt after 1 tick CALL IPOKE(IBPR,ICOUNT) CALL IPOKE(ICSR,IRATE) !start the clock 00 30 I e 1,10000 !do nothing loop while 30 CONTINUE !waiting for one tick 000C) 00 0 0 0 85 ITCK = IEND + l WRITE(7,*)'ITICK = ',ITCK ---5et zero (0) volts on all D/A outputs DO 40 I - 0,6,2 ! CALL IPOKE(IDAR+I,IZERO) !zero out all 4 channels 40 CONTINUE ! ---Get the INPUT A/D channels WRITE(7,500) !initial header 50 WRITE(7,510) !get INPUT channel READ(5,540,ERR=50) INCH IF ((INCH.GT.15).OR.(INCH.LT.0)) GOTO 50 7O WRITE(7,530) !get D/A output channel READ(5,540,ERR=70) IDACH IF ((IDACH.LT.0).OR.(IDACH.GT.3)) GOTO 70 IDAR = IDAR + IDACH * 2 ---Print out the options menu 80 WRITE(7,550) !print the menu options CALL PAUSE(SO,ITCK) !wait for printing to finish ---Disable the keyboard I - IPEEK(IK5R) .AND. "177477lclear keyboard interrupt CALL IPOKE(IKSR,I) ---Wait for keyboard input 90 CALL IPOKE(IPDR,IBELL) !ring the bell 100 I = IPEEK(IKSR) .AND. "200 !test bit 7 of IKSR IF (I .EQ. 0) GOTO 100 IDATA = IPEEK(IKDR) !a key has been pressed ---Check for a valid character IF ((IDATA.NE.'I').AND.(IDATA.NE.'G').AND.(IDATA.NE.'5').AND. + (IDATA.NE.'R'I.AND.(IDATA.NE.'Q'II GOTO 90 CALL IPOKE(IPDR,IDATA) !print the proper key ---Enable the keyboard I - IPEEK(IK5R) .OR. "100 !set bit 6 of IKSR CALL IPOKE(IK5R,I) WRITE(7,*)' ' ---Go to the proper place IF (IDATA .EQ. 'I') GOTO 110 IF (IDATA .EQ. 'G') GOTO 120 IF (IDATA .EQ. '5') GOTO 130 IF (IDATA .EQ. 'R') GOTO 140 IF (IDATA .EQ. 'Q') GOTO 150 —--INITIALIZE the zero reference position 110 CALL IPOSN 00000000000000000000 0 0 0 86 GOTO 80 !return to menu ---GAIN settings for PID control law 120 CALL CGAINS GOTO 80 !return to menu ---5TATU5 of Zero position 130 WRITE(7,560) !initial header WRITEI7,570) IZREF GOTO 80 !return to menu ---RUN the active control routine 140 IF (IZREF .NE. 0) GOTO 14S WRITE(7,600) !zero references GOTO 80 !not initialized 145 CALL CONTRL GOTO 80 !return to menu ---QUIT the program 150 WRITE(7,590) !exit message CALL PAUSE(60,ITCK) CALL EXIT END ------- - -— CGAINS -------------------------- SUBROUTINE CGAINS This subroutine allows the user to input the control gains for the PID control law, and the first order filter parameters. ---Definition of parameters F() F(1) = PID proportional gain (Kp) F(2) - PID derivative gain (Kd) F(3) - PID integral gain (Ki) F(d) = posn. filter parameter (alfa) F(S) - posn. filter parameter (1-alfa) F(6) - deriv. filter parameter (alfa) F(7) - deriv. filter parameter (1-alfa) F(8) 8 Input scale factor F(9) = Output scale factor F(10)= Maximum allowable error ---Declare Variables REAL G(10),F(10),SN COMMON /MISC2/ F,G,IZREF 0 87 ---Format Statements 100 FORMAT(/,/,2X,'INPUT the PID Control Gains',/, + 2X,‘ ') 110 FORMAT(/.2X,'1) Input PID Controller gains', /,2X,'2) Input Position Filter Parameters', /,2X,'3) Input Derivative Filter Parameters', /,2X,'4) Input INPUT/OUTPUT scale factors', /,2X,'5) Input maximum allowable error EMAX', /,2X,'6) Check Status of all gains', /,2x,'7) EXIT to main routine. .') 120 FORMAT(/,2X,'Your Selection: ',$) 130 FORMAT(Il) 140 FORMAT(/.2X,A14,' - ',S) 150 FORMAT(/.2X,'Enter I of Samples/Filter Time Const.: '.$) 160 FORMAT(/,2X.'PID Controller Gains:',/./. + 5X,'Kp - ',F10.4,/, + 5X,'Kd 8 ',F10.4,/. + 5X,'Ki = ',F10.4) 170 FORMAT(/,2X,A10,' Filter Parameters:',/. + 5X,'#5amples/Time Constant a ',F5.1,/, + 5x,' alpha = ',F10.4,/, + 5X,'1-alpha = ',F10.4) 180 FORMAT(/,2X,'Enter the ',A6,' scale [Volt/ ]: ',S) 190 FORMAT(/.2X,'Calibration Scale Factors:',/,/, + 5X,'Input Scale = ',F10.4,/, + 5X,'Output Scale - ',F10.4) 200 FORMAT(/,2X,'Enter the maximum allowable error: ',5) 210 FORMAT(/,2X,'Maximum Allowable Error:',/./. + SX,'Emax = ',F10.4) + + + + + + /. /. ---Main routine WRITE(7,100) 10 WRITE(7,110) !display menu options 20 WRITE(7,120) READ(5,130,ERR=20) IOP IF ((IOP.LT.1).OR.(IOP.GT.7)) GOTO 20 ---Goto the proper place GOTO(30,40,50,60,70,80,90) IOP ---Input new PID controller gains 30 WRITE(7,140) 'Kp [Volt/L] '!proportional gain READ(5,*,ERR=30) F(1) 32 WRITE(7,140) 'Kd [Volt-5/L]'!derivative gain READ(5,*,ERR=32) F(2) 34 WRITE(7,140) 'Ki [Volt/L-Sl'lintegral gain READ(5,*,ERR=34) F(3) GOTO 10 !return to CGAINS menu ---Input the filter parameters 40 WRITE(7,150) !get position filter READ(5,*,ERR=40) F(4) F(4) = 1.0 / (1.0 + F(4)) F(S) = 1.0 - F(4) 0 88 GOTO 10 ' !return to CGAINS menu ---Input the derivative filter parameters 50 WRITE(7,150) !get derivative filter READ(5,*,ERR=50) F(6) F(6) ' 1.0 / (1.0 + F(6)) F(7) - 1.0 - F(6) GOTO 10 !return to CGAINS menu ---Input the INPUT/OUTPUT scale factors 60 WRITE(7,180) 'INPUT' !get input scale factor READ(5,*,ERR=60) F(8) 65 WRITE(7,180) 'OUTPUT' !get output scale factor READ(5,*,ERR865) F(9) GOTO 10 !return to CGAINS menu ---Input the maximum allowable errror EMAX 7o WRITE(7,200) !get EMAX READ(S,*,ERR=70) F(10) GOTO 10 ---Check the status of all gains 80 WRITE(7,160) (F(I),I=1,3) SN - (1.0 / F(4II - 1.0 WRITE(7,170) 'Position ',SN,F(4),F(5) 5N - (1.0 / F(6)) - 1.0 WRITE(7,170) 'Derivative',5N,F(6),F(7) WRITE(7,190) F(8),F(9) WRITEI7,210) F(10) GOTO 10 !return to CGAINS menu ---EXIT to the main program 90 RETURN END 0 00000 SUBROUTINE CONTRL This is the control subroutine which sets up the clocked interrupt service control routine, and starts the control action. REAL G(lO),F(lO),FREQ,TICK,PER REAL EN,DEN,IEN,EO,DEO,IEO,VOLT 0 0 89 LOGICAL*1 IDATA.TIl(8),TI2(8) EXTERNAL RUNLSF !to use as subroutine arg. COMMON /ADDR1/ IKSR,IKDR,INCH,ITCK COMMON /ADDR2/ IPSR,IPDR,IBPR,IADVEC,IADERR COMMON /ADDR3/ ICSR,IADR,IDAR,IADBUF COMMON /MISC1/ ICHAN,PER,IFAST COMMON /MISC2/ F,G,IZREF COMMON /STATE/ EO,DEO,IEO,EN,DEN,IEN,VOLT ---Format Statements 100 FORMAT(/./.2X,'RUN the Clocked Control Routine', + /'2x'I _l) 110 FORMAT(/.2X,'Enter the SAMPLING FREQUENCY (Hz.): ',5) 120 FORMAT(/,2X,'ERROR. . .[Sampling frequency is too high.]') 130 FORMAT(/,2X,‘ERROR. . .[Sampling frequency is too low.]') 140 FORMAT(/.2X,'ERROR. . .[Sampling rate TOO HIGH for system.]') 150 FORMAT(/,2X,'Actual sampling freq. = ',E11.4,' Hz.') 160 FORMAT(/,2X,'Actual sample period - ',E11.4,' 5ec.',/,/) 170 FORMAT(/.2X,'Press ANY KEY to STOP active control. . .') 180 FORMAT(/,2X,'Are you ready to begin control [Y,N]. . . ',S) 190 FORMAT('+',5X,'Error: ',F10.5,5X,'dError: ',F10.5, + 5X,'Output: ',F10.5) 220 FORMAT('+',2X,'Turn OFF the Line Time Clock. . .') 230 FORMAT(/,2X,'ERROR. . .[System too far from origin]') 240 FORMAT(/.2X,'Is system near the origin [Y,N]. . . ',S) ---Initialize Variables IBELL = "7 IPR a 7 !highest priority interrupt IFAST = 0 !reset too-high flag FMAX - 200.0 !set max. sampling frequency VOLT = 0 !initial voltage ---Convert real gains to controller gains G(l) = - F(1) * F(9) / F(BIIKP .\ G(2) = - F(2) * F(9) / F(8)!Kd PID controller gains C(3) = - F(3) * F(9) / F(8)!Ki / G(4) = F(4) !position G(5) = F(5) !filter G(6) = F(6) !derivative G(7) = F(7) !filter C(10) = F(10) * F(8) * 204.8!EMAX C(10) = ABS(G(10)I ---5et up clocked ISR WRITE(7,100) !intial header 10 I = INTSET(IADVEC,IPR,1,RUNL5F)!attach to A/D done vector IF (I .EQ. 0) GOTO 20 WRITE(7,*) 'INTSET error -- A/D Vector, CODE = ',I GOTO 10 20 I = INTSET(IADERR,IPR,2,RUNL5F)!attach to A/D error vector 0 0 0 00 30 35 40 45 50 55 60 IF (I .EQ. 0) coro 30 WRITE(7,*) 'INTSET error -- A/D Error Vector, CODE = ',I GOTO 20 —--Input the sampling frequency WRITE(7,110) !get the sampling frequency READ(5,*,ERR=30) FREQ IF ((FREQ.LE.FMAX).AND.(FREQ.GT.0.)) GOTO 35 WRITE(7,120) !sampling frequency too high GOTO 30 ---Calculate the best base clock rate IR 8 1 !start at highest clock rate TICK = (10.0**(7-IR))/FREQ IF (TICK .LT. 3276?.) GOTO 45linteger out of range IR - IR + 1 !next lower base frequency IF (IR .LE. 7) GOTO 40 WRITE(7,130) !sampling frequency too low GOTO 30 ---Calculate the ticks for IBPR ITICK = IFIX(-l.0*(TICK+0.5))!nearest integer ---Calculate the actual sampling frequency and period FREQ - (10.0**(7-IR))/FLOAT(-ITICK) PER = 1.0 / FREQ WRITE(7.150) FREQ WRITE(7,160) PER ---5et up the clock status and A/D status registers IRATE = (IR * 8) + 3 ICHAN = "40140 + INCH * (2**8)!set up input sample -—-Make sure the line time clock is turned off CALL TIME(TIl) !get the first time CALL PAUSE(60,ITCK) !wait for 1 second CALL TIME(TIZ) !get the second time IF ((TIl(8).EQ.T12(8)).AND.(TI1(7).EQ.TIZ(7))) GOTO 60 CALL IPOKE(IPDR,IBELL) WRITE(7,220) !message to user GOTO 55 ---Ready to begin active control WRITE(7,170) !user instructions WRITE(7,180) !stop message CALL PAUSE(40,ITCK) !wait for terminal ---Disable keyboard and wait for input [Y.N] 0 0 65 70 72 74 80 91 I = IPEEK(IKSR) .AND. "177677lreset bit 6 of IKSR CALL IPOKE(IKSR,I) I = IPEEKIIKSR) .AND. "200 !test bit 7 of IKSR IF (I .EQ. 0) GOTO 65 IDATA = IPEEK(IKDR) la key has been pressed IF ((IDATA.EQ.'Y').OR.(IDATA.EQ.'N')) CALL IPOKE(IPDR,IDATA) IF (IDATA .EQ. 'Y') GOTO 70 IF (IDATA .EQ. 'N') GOTO 85 CALL IPOKE(IPDR,IBELL) GOTO 65 ---0isable the printer I = IPEEKlIPSR) .AND. "177677lreset bit 6 of IPSR CALL IPOKE(IPSR,I) ---5et up initial conditions for the system CALL SAMPLE(INCH,INVAL,IFLAG) !get initial position EN ‘ (IZREF - INVAL) !initial error DEN = 0.0 IEN = 0.0 ---Check for initial state too far from origin IF (ABS(EN) .LE. C(10)) GOTO 74 IDATA = IPEEK(IKDR) !reset the keyboard status I = IPEEK(IPSR) .OR. "100 !set bit 6 of IPSR CALL IPOKE(IPSR,I) WRITE(7,*)' ' WRITEI7.240) GOTO 65 ---5tart clock and ISR CALL IPOKE(IBPR,ITICK) CALL IPOKE(ICSR,IRATE) !start clock ticking CALL IPOKE(IADR,ICHAN) !start A/D conversion -—-Wait for keyboard input or too-fast error I = IPEEK(IPSR) .OR. "100 !set bit 6 of IPSR CALL IPOKE(IPSR,I) WRITE(7,*) ' ' WRITE(7,*) ' ' WRITE(7,*) ' ' WRITE(7,*) ' ' IF (IFAST .NE. 0) GOTO 85 WRITE(7,190) EN,DEN,VOLT I = IPEEK(IKSR) .AND. "200 !test bit 7 of IKSR IF (I .EQ. 0) GOTO 80 ---Enable the keyboard and printer 92 C 85 I = IPEEK(IK5R) .OR. "100 !set bit 6 of IKSR CALL IPOKE(IKSR,I) C I - IPEEK(IP5R) .OR. "100 !set bit 6 of IPSR CALL IPOKE(IPSR,I) C C ---Check for too-fast error C IF (IFAST .EQ. 0) GOTO 95 IF (IFAST .EQ. 2) GOTO 90 C WRITE(7,*)' ' WRITE(7,140) !print an error message GOTO 95 C 90 WRITE(7,*)' ' WRITE(7,230) !print an error message C C ---Turn off clock and A/D converter C 95 CALL IPOKE(IADR,O) CALL IPOKE(ICSR,0) WRITE(7,*)' ' C C ---Zero output voltage C IZERO a "4000 CALL IPOKE(IDAR,IZERO) C RETURN END C C _______________________________ -.. ------ c---— -------- RUNLSF -- ---- C SUBROUTINE RUNL5F(ID) C C This is an interrupt service routine to get clocked samples C of the input signal and compute the PID output control voltage. ---Definition of parameters G() 6(1) = proportional gain G(Z) = derivative gain C(3) = integral gain G(4) = filter parameter (alfa) G(S) = filter parameter (l-alfa) G(6) = deriv. filter parameter (alfa) G(7) - deriv. filter parameter (1-a1fa) 00000000000 REAL G(lO),F(10),VOLT,PER,VMAX REAL EO,DEO,IEO,EN,DEN,IEN COMMON /ADDR3/ ICSR,IADR,IDAR,IADBUF COMMON /MISC1/ ICHAN,PER,IFAST COMMON /MISC2/ F,G,IZREF COMMON ISTATE/ EO,DEO,IEO,EN,DEN,IEN,VOLT DATA VMAX/2047.5/ 0 --—Select entry: ID = 1 . . . A/D Sample Done [PID control] 0 00 00 0 00000 00 0000000 000 10 20 93 ID = 2 . . . A/D Sample Error GOTO (10,20) ID WRITE(7,*) 'RUNLSF entry error, 10 = ',ID CALL EXIT <<<<< Calculate PID Errors >>>>> -—-Get the new input sample EN = G(4)*(IZREF - IPEEK(IADBUF)) + G(5)*EO CALL IPOKE(IADR,ICHAN) lreset A/D status register ---Calcu1ate derivative of the input error DEN - G(6)*(EN - EO)/PER + G(7)*DEO -—-Ca1culate the integral of the input error IEN = IEO + EN * PER ---5hift new variables into the old variables E0 = EN DEO - DEN IEO = IEN —--Use these approximated errors for the PID control law <<<<< Calculate PID Control Law >>>>> ---Calculate the PID control law VOLT 8 G(1)*EN + G(2)*DEN + G(3)*IEN ---Check for voltage out of range IF (VOLT .LT. -VMAX) VOLT - -VMAX IF (VOLT .GT. VMAX) VOLT = VMAX IVOLT = IFIX(VOLT + 2048.0) CALL IPOKE(IDAR,IVOLT) !apply the voltage RETURN ---A/D Sample error CALL IPOKE(IADR,O) !turn off A/D CALL IPOKE(ICSR,0) !turn off RTC IFAST = 1 !set the too-fast flag IF (ABS(EN) .GT. C(10)) IFAST = 2 !out of range RETURN END --------------------------------- PID ---------------------------- 94 B.3 PD Control with Nonlinear Stick-Slip Compensation The PD control algorithm is implemented in the following program. A single input voltage is read and filtered with a first order digital filter. The derivative of this filtered signal is computed and filtered as well at user specified time constants. The nonlinear compensation force is also calculated. The compensated PD control signal is then computed and applied. This program allows the user to interactively tune all parameters. C C ——————————————————— ___ .. C --------------------- PD - C PROGRAM PD C C This is the main program which accesses all the control C subroutines for P0 control. C C Written by: STEVE C. SOUTHWARD C C ---Program Calling Tree C C PD C IPOSN C SAMPLE C CGAINS C CONTRL C RUNLSF C PAUSE C SETICK C C ---Declare Variables C REAL C(13).F(13) C LOGICAL*1 IDATA C DIMENSION IDLIST(21) C EXTERNAL SETICK C COMMON /ADDR1/ IKSR,IKDR,INCH,ITCK COMMON /ADDR2/ IPSR,IPDR,IBPR,IADVEC,IADERR COMMON /ADDR3/ ICSR,IADR,IDAR,IADBUF COMMON /MISC2/ F,G,IZREF COMMON /TIME/ I,IEND C C ---Initialize the zero reference C DATA IZREF/O/ C C ---Initial PD control gains, and other parameters C F(1) = 20.0 !Kp \ F(2) = 0.5 !Kd / PD controller gains 0 0 F(3) = 1.0/(1.0 + 0.0) F(4) = 1.0 - F(3) F(5) = 1.0/(1.0 + 100.0) F(6) = 1.0 - F(5) F(7) = 19.895 F(8) = 0.75936 F(9) 3 0.4 F(10) 2.0 F(11) - -2.0 F(12) ' 10.0 F(l3) ‘ 1.0 ---Set up address registers IPSW - "102 IKSR - ”177560 IKDR = IKSR + 2 IPSR - IKSR + 4 IPDR - IKSR + 6 ICSR - "170420 IBPR - ICSR + 2 ICKV ' "440 IDAR ' "170440 IADR = "170400 IADBUF ' IADR + 2 IADVEC - "400 IADERR - IADVEC + 4 IBELL - "7 IZERO 8 "4000 ---Protect System from crash IDLIST(1) = IKSR IDLIST(2) = IPEEK(IK5R) IDLIST(3) = IPSR IDLIST(4) = IPEEK(IP5R) IDLIST(S) - IDAR IDLIST(6) - IZERO IDLIST(7) = IDAR + 2 IDLIST(B) - IZERO IDLIST(9) - IDAR + 4 IDLIST(10) - IZERO IDLIST(11) a IDAR + 6 IDLIST(12) - IZERO IDLIST(13) = IADR IDLIST(14) = o IDLIST(15) = IADVEC IDLIST(16) = IPEEK(IADVEC) IDLIST(17) - ICSR IDLIST(18) - IPEEK(IC5R) IDLIST(19) - IPsw IDLIST(ZO) - IPEEK(IPSW) IDLIST(21) = 0 CALL DEVICE(IDLIST) ---Format Statements 95 !position !filter !derivative !filter !input scale factor !output scale factor !maximum allowable error !Fsp (+ve static friction) !an (-ve static friction) !Est (small +ve number) !Eps (small +ve number) !processor status word !keyboard status register !keyboard data register !printer status register !printer data register !real time clock status reg. !clock buffer preset reg. !clock vector address !D/A base output addr.: ch.0 !A/D control status register !A/D buffer register !A/D done interrupt vector !A/D error interrupt vector !bell tone !zero value for D/A converter !keyboard status !printer status !zero D/A channel 0 !zero D/A channel 1 !zero D/A channel 2 !zero D/A channel 3 !A/D status !A/D done interrupt vector !clock status !processor status !end of list marker FORMAT(/././.2X,'Discrete Time PD Control', + /,2X,' -------------------- _____ I) 0 0000 0 00 0 0 510 530 540 550 560 570 590 600 10 20 3O 4O 50 70 80 96 FORMAT(/.2X,'Enter the INPUT A/D channel: ',$) FORMAT(/.2X,'Enter the OUTPUT D/A channel: ',$) FORMAT(IZ) FORMAT(/,/,1OX,' ----- MAIN MENU ----- ', + /,/,2X,'(I)nitialize the zero reference position', + /,/,2X,'(G)ain settings for P0 control & filters', + /./,2X,'(5)tatus of Zero position', + /./.2X,'(R)un the active PD controller', + /./,2X.'(Q)uit the program', + /./.2X,'Choose one of the above. . . ',$) FORMAT(/./,2X,'5TATUS of Zero Reference. . .', + /.2X.' --- r') FORMAT(/,2X,'INPUT zero reference: ',IG) FORMAT(/./,2X,'Exiting Control Routine. . .',5(/)) FORMAT(/,/,2X,'Zero reference position NOT initialized') ---Calibrate the PAUSE subroutine IEND - 10000 !set initial value IPR a 7 !highest priority interrupt I = INTSET(ICKV,IPR,1,5ETICK) !attach to RTC vector IF (I .EQ. 0) GOTO 20 WRITE(7,*)'INT5ET error -- RTC vector, CODE = ',I GOTO 10 IRATE = "111 !set up clock to generate ICOUNT - —16667 !an interrupt after 1 tick CALL IPOKE(IBPR,ICOUNT) CALL IPOKE(ICSR,IRATE) !start the clock DO 30 I = 1,10000 !do nothing loop while CONTINUE !waiting for one tick ITCK = IEND + 1 WRITE(7,*)'ITICK = ',ITCK ---5et zero (0) volts on all D/A outputs DO 40 I - 0,6,2 ! CALL IPOKE(IDAR+I,IZERO) !zero out all 4 channels CONTINUE I ---Get the INPUT A/D channels WRITE(7,500) !initial header WRITE(7,510) !get INPUT channel READ(5,540,ERR=50) INCH IF ((INCH.GT.15).OR.(INCH.LT.O)) GOTO 50 WRITE(7,530) !get D/A output channel READ(5,540,ERR=70) IDACH IF ((IDACH.LT.0).OR.(IDACH.GT.3)) GOTO 70 IDAR = IDAR + IDACH * 2 ---Print out the options menu WRITE(7,550) !print the menu options CALL PAUSE(50,ITCK) !wait for printing to finish ---Disable the keyboard 0 0 0 0 97 I = IPEEK(IK5R) .AND. "177477lclear keyboard interrupt CALL IPOKE(IKSR,I) ---Wait for keyboard input 90 CALL IPOKE(IPDR,IBELL) !ring the bell 100 I a IPEEK(IKSR) .AND. "200 !test bit 7 of IKSR IF (I .EQ. 0) GOTO 100 IDATA - IPEEK(IKDR) !a key has been pressed ---Check for a valid character IF ((IDATA.NE.'I').AND.(IDATA.NE.'G').AND.(IDATA.NE.'S').AND. + (IDATA.NE.'R').AND.(IDATA.NE.'Q')) GOTO 90 CALL IPOKE(IPDR,IDATA) !print the proper key ---Enable the keyboard I = IPEEK(IKSR) .OR. "100 !set bit 6 of IKSR CALL IPOKE(IKSR,I) WRITE(7,*)' ' ---Go to the proper place IF (IDATA .EQ. 'I') GOTO 110 IF (IDATA .EQ. 'G') GOTO 120 IF (IDATA .EQ. '5') GOTO 130 IF (IDATA .EQ. 'R') GOTO 140 IF (IDATA .EQ. 'Q') GOTO 150 ---INITIALIZE the zero reference position 110 CALL IPOSN GOTO 80 !return to menu ---GAIN settings for P0 control law 120 CALL CGAINS GOTO 80 !return to menu ---5TATUS of Zero position 130 WRITE(7,560) !initial header WRITE(7,570) IZREF GOTO 80 !return to menu --—RUN the active control routine 140 IF (IZREF .NE. 0) GOTO 145 WRITE(7,600) !zero references GOTO 80 !not initialized 145 CALL CONTRL GOTO 80 !return to menu ---QUIT the program 00000000000000000000000 0 0 00 ------ - CGAINS -----------—- 98 150 WRITE(7,590) !exit message CALL PAUSE(60,ITCK) CALL EXIT END SUBROUTINE CGAINS This subroutine allows the user to input the control gains for the PD control law, and the first order filter parameters. ---Definition of parameters F(I F(1) = P0 proportional gain (Kp) F(2) - PD derivative gain (Kd) F(3) = posn. filter parameter (alfa) F(4) = posn. filter parameter (l-alfa) F(S) - deriv. filter parameter (alfa) F(6) = deriv. filter parameter (1-alfa) F(7) - Input scale factor F(8) - Output scale factor F(9) 2 Maximum allowable error F(10)- Fsp (+ve static friction) F(11)- an (-ve static friction) F(12)= Small positive number in control law (Eps) F(13)- Small positive number for error bound ---Declare Variables REAL C(13).F(13),SN COMMON /MISC2/ F,G,IZREF ---Format Statements 100 FORMAT(/,/,2X,'INPUT the PD Control Gains',/, + 2x.' - ') 110 FORMAT(/,2X,'1) Input P0 Controller gains', /,2X,'2) Input Position Filter Parameters', /,2X,'3) Input Derivative Filter Parameters', /,2X,'4) Input INPUT/OUTPUT scale factors', /,2X,'5) Input maximum allowable error EMAX', /,2X,'6) Input Stiction compensation gains', /,2X,'7) Check Status of all gains', /,2X,'8) EXIT to main routine. .') 120 FORMAT(/,2X,'Your Selection: ',$) 130 FORMAT(Il) 140 FORMAT(/,2X,A14,' - ',$) 150 FORMAT(/.2X,'Enter I of Samples/Filter Time Const.: ',$) 160 FORMAT(/,2X,'PD Controller Gains:',/,/, + 5X,'Kp - ',F10.4,' [N/m]',/. + 5X,'Kd - ',F10.4,' [N-s/ml') 170 FORMAT(/.2X,A10,' Filter Parameters:',/,/, + 5x,'ISamples/Time Constant = ',F5.1,/ / + 5X,‘ alpha - ',F10.4,/, + 5X,'1-alpha - ',F10.4) 180 FORMAT(/,2X,'Enter the ',A6,' scale: ',$) + + + + + + + I I 0 0 0 99 190 FORMAT(/,2X,'Calibration Scale Factors:',/,/. + 5X,'Input Scale = ',F10.4,' [Volt/m]',/, + 5X,'Output Scale - ',F10.4,' [N/Voltl') 200 FORMAT(/,2X,'Enter the maximum allowable error [m]: ',$) 210 FORMAT(/.2X,'Maximum Allowable Error:'././, + 5X,'Emax - ',F10.4,' [m]') 220 FORMAT(/,2X,'Stiction Compensation Gains:',/,/. 5X,'Fsp - ',F10.4,' [N]',/, 5X,'an = ',F10.4,' [N]',/. 5X,'Est - ',F10.4,' [bits]',/, 5X,'EPS - ',F10.4,' [bits]') + + + + ---Main routine WRITE(7,100) 10 WRITE(7,110) !display menu options 20 WRITE(7,120) READ(5,130,ERR=20) IOP IF ((IOP.LT.1).OR.(IOP.GT.8)) GOTO 20 ---Goto the proper place GOTO(30,40,50,60,70,80,85,90) IOP ---Input new PD controller gains 30 WRITE(7,140) 'Kp [N/m] ' !proportional gain READ(5,*,ERR=30) F(1) 32 WRITE(7,140) 'Kd [N‘s/m1' !derivative gain READ(5,*,ERR=32) F(2) GOTO 10 !return to CGAINS menu ---Input the filter parameters 40 WRITE(7,150) !get position filter READ(5,*,ERR=40) F(3) F(3) = 1.0 / (1.0 + F(3)) F(4) = 1.0 - F(3) GOTO 10 !return to CGAINS menu ---Input the derivative filter parameters 50 WRITE(7,150) !get derivative filter READ(5,*,ERR=50) F(5) F(5) = 1.0 / (1.0 + F(5)) F(6) = 1.0 - F(5) GOTO 10 !return to CGAINS menu ---Input the INPUT/OUTPUT scale factors 60 WRITE(7,180) 'Input Cx [V/m]'!get input scale factor READ(5,*,ERR=60) F(7) 65 WRITE(7,180) 'Outpt Cu [N/VJ'lget output scale factor READ(5,*,ERR=65) F(8) GOTO 10 !return to CGAINS menu 0 0 0 00000 100 ---Input the maximum allowable errror EMAX 70 WRITE(7,200) !get EMAX READ(5,*,ERR=70) F(9) GOTO 10 ---Input the stiction compensation gains 80 WRITE(7,140) 'Fsp [Nl' !get Fsp READ(5,*,ERR=80) F(10) 81 WRITE(7,140) 'an [Nl' !get an READ(5,*,ERR=81) F(ll) 82 WRITE(7,140) 'Est [bits]' !get Epsilon for X READ(5.*,ERR=82) F(12) 83 WRITE(7,140) 'EPS [bits]' !get precision bound READ(5,*,ERR=83) F(13) GOTO 10 ---Check the status of all gains 85 WRITE(7,160) (F(I),I=1,2) SN ' (1.0 / F(3)) - 1.0 WRITE(7,170) 'Position ',SN,F(3),F(4) 5N .. (1.0 / F(5)) - 1.0 warrst7,170) 'Derivative',SN,F(5),F(6)I WRITE(7,190) F(7),F(8) WRITE(7,210) F(9) WRITE(7,220) F(10),F(11),F(12),F(13) GOTO 10 !return to CGAINS menu ---EXIT to the main program 90 RETURN END SUBROUTINE CONTRL This is the control subroutine which sets up the clocked interrupt service control routine, and starts the control action. REAL C(13),F(l3),FREQ,TICK,PER REAL EN,DEN,EO,DEO,VOLT LOGICAL*1 IDATA,TIl(8),TI2(8) EXTERNAL RUNLSF !to use as subroutine arg. 0 0 100 110 120 130 140 150 160 170 180 185 186 187 190 220 230 240 COMMON COMMON COMMON COMMON COMMON COMMON 101 /ADDR1/ IKSR,IKDR,INCH,ITCK /ADDR2/ IPSR,IPDR,IBPR,IADVEC,IADERR IADDR3/ ICSR,IADR,IDAR,IADBUF /MISC1/ ICHAN,PER,IFAST /MISC2/ F,G,IZREF [STATE/ EO,DEO,EN,DEN,VOLT -—-Format Statements FORMAT(/,/,2X,'RUN the Clocked Control Routine', + /,2x, FORMAT(/.ZXt FORMAT(/.ZX. FORMAT(/IZX. FORMAT(/.ZXo FORMAT(/,ZX, FORMAT(/,ZX, FORMAT(/.ZX, FORMAT(/.ZX: FORMAT(/.ZX. /.2X, FORMAT(/,ZX, + FORMAT( Il) ') 'Enter the SAMPLING FREQUENCY (Hz.): ',$) 'ERROR. . .[Sampling frequency is too high.]') 'ERROR. . .[Sampling frequency is too low.]') 'ERROR. . .[Sampling rate TOO HIGH for system.]') 'Actual sampling freq. - ',E11.4,' Hz.') 'Actual sample period - ',E11.4,' 5ec.',/,/) 'Press ANY KEY to STOP active control. . .') 'Are you ready to begin control [Y,N]. . . ',$) '1) PD control only', '2) PD w/Stiction compensation') 'Choose one [1,2]: ',$) FORMAT('+',5X,'Error: ',F10.4,5X,'Output: ',F10.4) FORMAT('+',2X,'Turn OFF the Line Time Clock. . .') FORMAT(/,2X,'ERROR. . .[System too far from origin]') FORMAT(/,2X,'Is system near the origin [Y,N]. . . ',$) ---Initialize Variables IBELL = IPR - 7 IFAST 8 FMAX = 250.0 VOLT 8 "7 0 0 !highest priority interrupt !reset too-high flag !set max. sampling frequency !initial voltage ---Convert real gains to controller gains G(l) - F(1) / (F(8) * F(7)) !KP \ G(2) = F(2) / (F(8) * F(7)) !Kd / PD controller gains C(3) = F(3) !position G(4) = F(4) !filter G(S) = F(5) !derivative G(6) = F(6) !filter G(9) = F(9) * F(7) * 204.8 !EMAX G(9) - ABS(G(9)) G(lO) = - (F(10) / F(1)) * F(7) * 204.8 C(11) = - (F(11) / F(1)) * F(7) * 204.8 C(12) = ABS(F(12)) C(10) = G(lO) - C(12) C(11) = C(11) + C(12) C(13) - ABS(F(13)) ---5et up clocked ISR WRITE(7.100) WRITE(7.185) !intial header 0 0 5 10 15 20 30 35 40 45 50 102 WRITE(7,186) !choose desired control alg. READ(5,187,ERR=5) IOP IF ((IOP.NE.1).AND.(IOP.NE.2)) GOTO 5 GOTO (10,15) IOP I - INTSET(IADVEC,IPR,1,RUNL5F)!attach to A/D done vector IF (I .EQ. O) GOTO 20 WRITE(7,*) 'INTSET error -- A/D Vector, CODE = ',I GOTO 10 I = INTSET(IADVEC,IPR,2,RUNLSF)!attach to A/D done vector IF (I .EQ. 0) GOTO 20 WRITE(7,*) 'INTSET error -- A/D Vector, CODE = ',I GOTO 15 I - INTSET(IADERR,IPR,3,RUNL5F)!attach to A/D error vector IF (I .EQ. 0) GOTO 30 WRITE(7,*) 'INTSET error -- A/D Error Vector, CODE = ',I COTO 20 ---Input the sampling frequency WRITE(7,110) !get the sampling frequency READ(5,*,ERR=30) FREQ IF ((FREQ.LE.FMAX).ANO.(FREQ.GT.0.)) GOTO 35 WRITE(7,120) !sampling frequency too high GOTO 30 ---Calculate the best base clock rate IR = 1 !start at highest clock rate TICK = (10.0**(7-IR))/FREQ IF (TICK .LT. 32767.) GOTO 45!integer out of range IR = IR + 1 !next lower base frequency IF (IR .LE. 7) GOTO 40 WRITE(7,130) !sampling frequency too low GOTO 30 -‘-Calcu1ate the ticks for IBPR ITICK = IFIX(-1.0*(TICK+O.5))!nearest integer -—-Calculate the actual sampling frequency and period FREQ ' (10.0**(7-IR))/FLOAT(-ITICK) PER = 1.0 / FREQ G(5) = G(5) / PER WRITE(7,150) FREQ WRITE(7,160) PER ---Set up the clock status and A/D status registers IRATE = (IR * 8) + 3 ICHAN = "40140 + INCH * (2**8)!set up input sample 0 0 0 0 55 60 65 70 72 74 103 ---Make sure the line time clock is turned off CALL TIME(T11) !get the first time CALL PAU5E(60,ITCK) !wait for 1 second CALL TIME(TIZ) !get the second time IF ((TIl(8).EQ.T12(8)).AND.(TIl(7).EQ.TIZ(7))) GOTO 60 CALL IPOKE(IPDR,IBELL) WRITE(7,220) !message to user GOTO 55 ---Ready to begin active control WRITE(7,170) !user instructions WRITE(7,180) !stop message CALL PAUSE(40,ITCK) !wait for terminal -—-Disab1e keyboard and wait for input [Y.NI I = IPEEK(IKSR) .AND. "177677 !reset bit 6 of IKSR CALL IPOKE(IKSR,I) I = IPEEK(IKSR) .AND. "200 !test bit 7 of IKSR IF (I .EQ. 0) GOTO 65 IDATA 2 IPEEK(IKDR) la key has been pressed IF ((IDATA.EQ.'Y').OR.(IDATA.EQ.'N')) CALL IPOKE(IPDR,IDATA) IF (IDATA .EQ. 'Y') GOTO 70 IF (IDATA .EQ. 'N') COTO 85 CALL IPOKE(IPDR.IBELL) GOTO 65 ---Disable the printer I = IPEEK(IPSR) .AND. "177677lreset bit 6 of IPSR CALL IPOKE(IPSR,I) ---5et up initial conditions for the system CALL SAMPLE(1NCH,INVAL,IFLAG) !get initial position EN = (IZREF - INVAL) linitial error DEN = 0.0 ---Check for initial state too far from origin IF (ABS(EN) .LE. G(9)) COTO 74 IDATA = IPEEK(IKDR) !reset the keyboard status I 8 IPEEK(IPSR) .OR. "100 !set bit 6 of IPSR CALL IPOKE(IPSR,I) WRITE(7,*)' ' WRITE(7,240) GOTO 65 ---Start clock and ISR CALL IPOKE(IBPR.ITICKI 104 CALL IPOKE(ICSR,IRATE) !start clock ticking CALL IPOKE(IADR,ICHAN) !start A/D conversion C C --—Wait for keyboard input or too-fast error C I = IPEEK(IP5R) .OR. "100 !set bit 6 of IPSR CALL IPOKE(IPSR,I) C WRITE(7,*) ' ' WRITE(7,*) ' ' WRITE(7,*) ' ' WRITE(7,*) ' ' 80 IF (IFAST .NE. 0) GOTO 85 WRITE(7,190) EN,VOLT C I = IPEEK(IK5R) .AND. "200 !test bit 7 of IKSR IF (I .EQ. 0) GOTO 80 C C --—Enable the keyboard and printer C 85 I = IPEEK(IKSR) .OR. "100 !set bit 6 of IKSR CALL IPOKE(IKSR,II C I - IPEEK(IPSR) .OR. "100 !set bit 6 of IPSR CALL IPOKE(IPSR,I) C C ---Check for too-fast error C IF (IFAST .EQ. 0) GOTO 95 IF (IFAST .EQ. 2) GOTO 90 C WRITE(7,*)' ' WRITE(7,140) !print an error message GOTO 95 C 90 WRITE(7,*)' ' WRITE(7,230) !print an error message C C ---Turn off clock and A/D converter C 95 CALL IPOKE(IADR,O) CALL IPOKE(ICSR,0) WRITE(7,*)' ' C C ---Zero output voltage C IZERO = "4000 CALL IPOKE(IDAR,IZERO) C RETURN END C C —————————————————————————————————————————————————————————————————— C- - --- RUNLSF -------------------------- C SUBROUTINE RUNL5F(ID) C C This is an interrupt service routine to get clocked samples C of the input signal and compute the PD output control voltage. C C ---Definition of parameters G() C 0000000 0 00000 00000 0 0 0 00 0000000 0 0 105 G(l) - proportional gain G(2) - derivative gain G(3) - filter parameter (alfa) G(4) filter parameter (1-a1fa) G(5) deriv. filter parameter (alfa) G(6) = deriv. filter parameter (l-alfa) REAL G(l3),F(l3),VOLT,PER,VMAX REAL EO,DEO,EN,DEN COMMON /ADDR3/ ICSR,IADR,IDAR,IADBUF COMMON /MISC1/ ICHAN,PER,IFAST COMMON /MISC2/ F,G,IZREF COMMON /STATE/ EO,DEO,EN,DEN,VOLT DATA VMAX/2047.5/ ---5elect entry: ID - 1 . . . A/D Sample Done [PD control] ID = 2 . . . A/D Sample Done [PD w/Comp] ID = 3 . . . A/D Sample Error IF (ID .EQ. 3) GOTO 90 <<<<< Calculate PD Errors >>>>> ---Get the new input sample EN = C(3)*(IZREF - IPEEK(IADBUF)) + G(4)*EO IF (ABS(EN) .GT. G(9)) GOTO 90!out of range CALL IPOKE(IADR,ICHAN) !reset A/D status register ---Calculate derivative of the input error DEN a G(5)*(EN - E0) + G(6)*DEO ---Shift new variables into the old variables E0 = EN DEO = DEN ---Use these approximated errors for the PD control law <<<<< Calculate PD Control Law >>>>> ---Calcu1ate the PD control law IF (ID .EQ. 2) GOTO 10 VOLT 8 G(l)*EN + G(2)*DEN GOTO 80 ---Compute the compensation force 10 IF ((‘EN.GE.G(10)).AND.(-EN.LE.G(11))I GOTO 20 VOLT = G(l)*EN + G(2)*DEN GOTO 50 20 IF (‘EN .LE. G(l3)) GOTO 30 VOLT = - G(l) * G(ll) + G(2)*DEN GOTO 50 30 IF (-EN .GE. -G(13)) GOTO 40 VOLT - - G(l) * C(10) + G(2)*DEN GOTO 50 0 106 40 CONTINUE VOLT = G(2)*DEN 50 CONTINUE ---Check for voltage out of range 80 IF (VOLT .LT. -VMAX) VOLT = -VMAX IF (VOLT .GT. VMAX) VOLT 8 VMAX IVOLT ' IFIX(VOLT + 2048.0) CALL IPOKE(IDAR,IVOLT) !apply the voltage RETURN ---A/D Sample error 90 CALL IPOKE(IADR,0) !turn off A/D CALL IPOKE(ICSR,0) !turn off RTC IFAST = 1 !set the too-fast flag IF (ABS(EN) .GT. G(9)) IFAST RETURN END B.4 Common Control Modules 2!out of range The following subroutines are common to both the PID and PD control programs. These are IPOSN, PAUSE, SETICK, and SAMPLE. 0000 0 0 C 0 C 0 ............. -- --- IPOSN SUBROUTINE IPOSN This subroutine allows the user to define the initial zero reference position of the input signal. REAL G(IO).F(10) LOGICAL*1 IDATA COMMON /ADDR1/ IKSR,IKDR,INCH,ITCK COMMON /MISC2/ F,G,IZREF —~-Format Statements 100 FORMAT('+',5X,'INPUT Sample: ',I4) 110 FORMAT(/./,2X,'INITIALIZE the Reference Position',/, 2X,‘ + + + + ---Main Routine 2X,'Move the system to its reference',/,/, 2X,'position to define the ZERO.'././. 2X,'Press "X" to define the INPUT reference -'./././) WRITE(7,110) !print the instructions 0 0 00000 0 0 0000 00 107 CALL PAUSE(100,ITCK) !wait for CRT ---Turn off the keyboard I = IPEEK(IKSR) .AND. "177477Eclear keyboard interrupt CALL IPOKE(IKSR,I) ---Start Sampling 10 CALL SAMPLE(INCH,INVAL,IFLAG)!get INPUT A/D sample IF (IFLAG .NE. 0) GOTO 20 IZREF = INVAL WRITE(7,100) IZREF CALL PAUSE(7,ITCK) !wait for 7 ticks ---Check for keyboard input I = IPEEK(IKSR) .AND. "200 !test bit 7 of IKSR IF (I .EQ. 0) GOTO 10 IDATA 8 IPEEK(IKDR) !a key has been pressed IF (IDATA .EQ. 'X') GOTO 20 GOTO 10 --—Turn on the keyboard 20 I - IPEEK(IKSR) .OR. "100 CALL IPOKE(IKSR,I) RETURN END SUBROUTINE PAUSE(NTICKS,ITCK) This subroutine is approximately calibrated to pause in ITCK intervals, as specified on input. There are 60 ticks in one second. DO 20 I = l,NTICKS DO 10 J = 1,ITCK !one tick per outer 10 CONTINUE !do loop step on (I) 20 CONTINUE RETURN END - SETICK ---- SUBROUTINE SETICK(ID) This is an interrupt service routine that is called 108 C when the real time programmable clock has generated C an interrupt after one tick (1/60 sec.). C COMMON /TIME/ I,IEND C IF (ID .EQ. 1) GOTO 10 !check for entry error C WRITE(7,*) 'SETICK entry error, ID = ',ID CALL EXIT C 10 IEND - I I = 10000 !reset do loop counter C RETURN END '0 —————————————— — —— -—-— ;- - -— SAMPLE ------------ .TITLE SAMPLE ; This is a FORTRAN compatible subroutine called with ; CALL SAMPLE(ICHAN,IVAL,IFLG) ; This subroutine gets a sample from the specified A/D converter ; channel ICHAN. The data read from the A/D is passed back in 'IVAL'. ; The flag IFLG is returned to the main program with the status: ; / 0 = NO ERRORS DETECTED ; IFLG = - 1 = CHANNEL NUMBER OUT OF RANGE (0-15) 3 \ 2 = A/D SAMPLING ERROR .GLOBL SAMPLE ;This makes SAMPLE a global symbol ‘0 LC=. ;Find current assembly location .=1000+LC ;Start assembling at location 1000 CSR=170400 ;Define A/D status register address DBR=C5R+2 ;Define A/D buffer register address SAMPLE: R5 points to the first address of a table of arguments @(R5) is the first data (= I of arguments passed) @2(R5) is the first argument passed ‘0 Q. ‘. MOV @2(R5),ICHAN ;Put first argument in ICHAN MOV @2(R5),TEMP ;Put first argument in TEMP also BIC #177760,ICHAN ;Clear bits 4 to 15 CLR FLG ;Clear the FLG CMP TEMP,ICHAN ;Check for ICHAN > 15 ENE ERRl ;Branch to ERRl if ICHAN > 15 START: SWAB ICHAN ;Swap high byte for low byte MOV ICHAN,@#C5R ;Set the A/D Control Status Register BIS #1,@#CSR ;Start the A/D conversion WAIT: BIT #200,@#CSR ;Test bit 7 for A/D done BEQ WAIT :Wait if not set MOV @IDBR,@4(R5) ;Put sampled data into argument table BIT #100000,@#C5R ;Test bit 15 for A/D sampling error BNE ERR2 ;Branch to ERR2 if bit 15 is set DONE: MOV FLG,@6(R5) ;Put flag back into argument table RTS PC I I ; ERROR HANDLING ROUTINES ERRl: BIS #1,FLG BR START ERR2: BIS #2,FLG BR DONE ; LOCAL VARIABLES ICHAN: .WORD 0 TEMP: FLG: I .WORD 0 .WORD 0 .END 109 ;Return from subroutine to calling program ;Set bit 0 of FLG ;Go back to start ;Set bit 1 of FLG ;Go back to DONE References References Amontons, G., 1699, Histoire de I’Academie Royale des Sciences avec les Memoires de Mathematique et de Physique, p. 206. Armstrong, B., 1988, “Friction: Experimental Determination, Modeling and Compensation,” IEEE International Conference on Robotics and Automation, Vol. 3. PP. 1422-1427. Barbashin, E.A., 1970, Introduction to the Theory of Stability, Wolters-Noordhoff Publishing, The Netherlands, pp. 66-68. Bowden, PP, and Tabor, D., 1950, The Friction and Lubrication of Solids, Oxford, Clarendon Press. Bowden, RR, and Tabor, D., 1964, The Friction and Lubrication of Solids Part 11, Oxford, Clarendon Press. Canudas, C., Astrom, KL, and Braun, K., 1987, “Adaptive Friction Compensation in DC—Motor Drives,” IEEE Journal of Robotics and Automation, Vol. RA-3, No. 6, pp. 681-685. Cheok, K.C., Hu, H., and Loh, N.K., 1988, “Modeling and Identification of a Class of Servomechanism Systems With Stick-Slip Friction,” Transactions of the ASME, Vol. 110. PP. 324-328. Coulomb, C.A., 1785, Memoires de Mathematique et de Physique de l’Academie Royale des Sciences, p. 161. Craig, J., 1988, Adaptive Control of Mechanical Manipulators, Addison-Wesley Publishing Company Inc. Davison, C. St. C., 1957, Wear, Vol. 1, pg. 155. Euler, L., 1750, Histoire de I’Academie Royale des Sciences at Belles-Lettres, 1748. Berlin, 1750, p. 122. 110 111 Filippov, A.F., 1964, “Differential Equations with Discontinuous Right-Hand Side,” American Mathematical Society Translations, Vol. 42, Series 2, pp. 199-231. Gilbart, J.W., and Winston, G.C., 1974, “Adaptive Compensation for an Optical Tracking Telesc0pe,”Automatica, Vol. 10, pp. 125-131. Gogoussis, A., and Donath, M., 1987, “Coulomb Friction Joint and Drive Effects in Robot Mechanisms,” IEEE, pp. 828-836. Hahn, W., 1963, Theory and Application of Lyapunov's Direct Method, Prentice-Hall, Englewood Cliffs, New Jersey, pg. 35. Hahn, W., 1967, Stability of Motion, Springer-Verlag, New York Inc., pp. 369-371. Halliday, D., and Resnick, R., 1978, Physics, Parts I & 2, Third Ed., John Wiley & Sons, New York, pp. 97-101. Hartog, J.P. Den, 1930, “Forced Vibrations with Combined Coulomb and Viscous Damping,” Transactions of the ASME, Vol. 53, pp. 107-115. Karnopp, D., 1985, “Computer Simulation of Stick-Slip Friction in Mechanical Dynamic Systems,” ASME Journal of Dynamic Systems, Measurement, and Control, Vol. 107. PP. 100-103. Kolston, P.J., 1988, “Modeling Mechanical Stick-Slip Friction Using Electrical Circuit Analysis,” Transactions of the ASME, Vol. 110, pp. 440-443. Kubo, T., Anwar, G., and Tomizuka, M., 1986, “Application of Nonlinear Friction Compensation to Robot Arm Control,” IEEE International Conference on Robotics and Automation, Vol. 2, pp. 722-727. McShane, EL, 1947, Integration, Princeton University Press, Princeton, pp. 188-217. Miller, E, 1977, College Physics, 4th Ed., Harcourt Brace Jovanovich, Inc., New York, pg. 88. Rabinowicz, B., 1951, “The Nature of the Static and Kinetic Coefficients of Friction,” Journal of Applied Physics, Vol. 22, No. 11, November, pp. 1373-1379. 112 Rouche, N., Habets, P., and Laloy, M., 1977, Stability Theory by Lyapunov’s Direct Method, Springer-Verlag, New York, pp. 345-346. Shaw, 8., 1986, “On the Dynamic Response of a System with Dry Friction,” Journal of Sound and Vibration, Vol. 108(2), pp. 305-325. Siljak, DD, 1969, Nonlinear Systems, The Parameter Analysis and Design, John Wiley & Sons, Inc., New York, pp. 539-544. Solncev, Yu. K., 1951, “On the Stability in the Sense of Lyapunov of the Equilibrium Positions of a System of two Differential Equations in the Case of Discontinuous Right-Hand Sides,” Moskovski Gos. Universitet, Uchenve Zapiski, Vol. 148, Mat. 4, pp. 144-180. Southward, S., 1986, A Modified Luenberger Observer for the Inverted Pendulum, MS Thesis, Michigan State University, East Lansing, Michigan. Southward, S., 1989, DIFFEQ User’s Guide, Double Precision, VMS Version 2.0. Tomlinson, G.R., and Hibbert, J .H., 1979, “Identification of the Dynamic Characteristics of a Structure with Coulomb Friction,” Journal of Sound and Vibration, Vol. 64, No. 2. PP. 233-242. Tou, J ., and Schultheiss, RM, 1953, “Static and Sliding Friction in Feedback Systems,” Journal of Applied Physics, Vol. 24, No. 9, pp. 1210-1217. Townsend, W.T., and Salisbury, J.K., 1987, “The Effect of Coulomb Friction and Stiction on Force Control,” IEEE International Conference on Robotics and Automation, Vol. 2, pp. 883-889. Tustin, A., 1947, “The Effects of Backlash and of Speed-Dependent Friction on the Stability of Closed-Cycle Control Systems,” Journal of the Institute of Electrical Engineers, Vol. 94(2A), pp. 143-151. Utkin, VI, 1977, “Variable Structure Systems with Sliding Modes,” IEEE Transactions on Automatic Control, Vol. AC-22, No. 2, April, pp. 212-222. Vidyasagar, M., 1978, Nonlinear Systems Analysis, Prentice-Hall, Englewood Cliffs, New Jersey. 113 Vinci, L. da, circa 1452-1519, S.K.M. III, The Forster Manuscripts, Library of the Victoria and Albert Museum, South Kensington, 1923-1930, 72 r. Vinci, L. da, circa 1452-1519, S.K.M. II, The Forster Manuscripts, Library of the Victoria and Albert Museum, South Kensington, 1923-1930, 133 r. Vinci, L. da, circa 1452-1519, Codex Atlanticus, Ambrosiana Library, Milan, 1894-1904, 198 v. Walrath, CD, 1984, “Adaptive Bearing Friction Compensation Based on Recent Knowledge of Dynamic Friction,” Automatica, Vol. 20, No. 6, pp. 717-727. Yang, S., and Tomizuka, M., 1987, “Adaptive Pulse Width Control for Precise Positioning Under Influence of Stiction and Coulomb Friction,” IEEE, pp. 188- 193. Yoshizawa, T., 1966, Stability Theory by Lyapunov's Second Method, The Mathematical Society of Japan, Gakujutsutosho Printing Co., Ltd. Tokyo, pp. 1-10. \Ftc 7”??? Ft pm B8 u. .N h_\€o¢+ MS coo I X I wloo HICH GRN STATE N V. l r u I LIBRARIES l I II “1MIUIWIIWIIWHI 0 7163 312930 572