Title Page
 Table of Contents
 List of Figures
 Review of existing results
 Development of standard model
 Existence conditions for limit...
 Limit cycle symmetry and stabi...
 Conclusions and recommendations...
 Biographical sketch

Title: Limit cycles for systems with nonlinear friction /
Full Citation
Permanent Link: http://ufdc.ufl.edu/UF00097391/00001
 Material Information
Title: Limit cycles for systems with nonlinear friction /
Physical Description: vii, 140 leaves : ill. ; 28 cm.
Language: English
Creator: James, Michael R
Publication Date: 1989
Copyright Date: 1989
Subject: Limit cycles   ( lcsh )
Servomechanisms   ( lcsh )
Nonlinear mechanics   ( lcsh )
Piecewise linear topology   ( lcsh )
Electrical Engineering thesis Ph. D
Dissertations, Academic -- Electrical Engineering -- UF
Genre: bibliography   ( marcgt )
non-fiction   ( marcgt )
Thesis: Thesis (Ph. D.)--University of Florida, 1989.
Bibliography: Includes bibliographical references (leaves 137-138)
Additional Physical Form: Also available on World Wide Web
Statement of Responsibility: by Michael R. James.
General Note: Typescript.
General Note: Vita.
 Record Information
Bibliographic ID: UF00097391
Volume ID: VID00001
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: alephbibnum - 001530507
oclc - 22331252
notis - AHE3883


This item has the following downloads:

limitcyclesforsy00jamerich ( PDF )

Table of Contents
    Title Page
        Page i
        Page ii
    Table of Contents
        Page iii
        Page iv
    List of Figures
        Page v
        Page vi
        Page vii
        Page 1
        Page 2
        Page 3
        Page 4
        Page 5
        Page 6
        Page 7
        Page 8
        Page 9
        Page 10
        Page 11
        Page 12
        Page 13
        Page 14
        Page 15
        Page 16
        Page 17
        Page 18
        Page 19
        Page 20
        Page 21
        Page 22
    Review of existing results
        Page 23
        Page 24
        Page 25
        Page 26
        Page 27
        Page 28
        Page 29
    Development of standard model
        Page 30
        Page 31
        Page 32
        Page 33
        Page 34
        Page 35
        Page 36
        Page 37
        Page 38
        Page 39
        Page 40
        Page 41
        Page 42
        Page 43
        Page 44
        Page 45
        Page 46
        Page 47
        Page 48
    Existence conditions for limit cycles
        Page 49
        Page 50
        Page 51
        Page 52
        Page 53
        Page 54
        Page 55
        Page 56
        Page 57
        Page 58
        Page 59
        Page 60
        Page 61
        Page 62
        Page 63
        Page 64
        Page 65
        Page 66
        Page 67
        Page 68
        Page 69
        Page 70
    Limit cycle symmetry and stability
        Page 71
        Page 72
        Page 73
        Page 74
        Page 75
        Page 76
        Page 77
        Page 78
        Page 79
        Page 80
        Page 81
        Page 82
        Page 83
        Page 84
        Page 85
        Page 86
        Page 87
        Page 88
        Page 89
        Page 90
        Page 91
        Page 92
        Page 93
        Page 94
        Page 95
        Page 96
        Page 97
        Page 98
        Page 99
        Page 100
        Page 101
        Page 102
        Page 103
        Page 104
        Page 105
        Page 106
        Page 107
        Page 108
        Page 109
        Page 110
        Page 111
        Page 112
        Page 113
        Page 114
        Page 115
    Conclusions and recommendations for further research
        Page 116
        Page 117
        Page 118
        Page 119
        Page 120
        Page 121
        Page 122
        Page 123
        Page 124
        Page 125
        Page 126
        Page 127
        Page 128
        Page 129
        Page 130
        Page 131
        Page 132
        Page 133
        Page 134
        Page 135
        Page 136
        Page 137
        Page 138
    Biographical sketch
        Page 139
        Page 140
        Page 141
        Page 142
        Page 143
Full Text








I am indebted to many people for their cooperation and

assistance over the long period of my doctoral studies.

First, my adviser, Dr. Thomas Bullock, and the other

members of my graduate committee, for their suggestions,

help and cooperation. Dr. V. M. Popov provided a very

helpful discussion on symmetry. Dean Eugene Chenette and

Art Zirger were very helpful during my coursework and

residence, as were many other UF personnel too numerous to


I received help from many fellow employees of Harris

Corporation: Pete Pitard, for his hard work in obtaining

financial assistance, Lee Almond, for his support, Kevin

Arter, for his review and suggestions, and the

Training/FEEDS and Human Resources personnel. Special

thanks go to Rich and Donna Phelan, whose enthusiastic

support solved many problems during the first years.

Most of all, I want to thank my wife, parents, and family

for their understanding and support during these six long

years. Without their sacrifices I could not have succeeded,

nor would it be meaningful without them.

All of these people share in this achievement, and I

deeply thank them for their support.


ACKNOWLEDGMENTS............................ .. ....... ii

LIST OF FIGURES ................................. ... v

ABSTRACT............................................ vi


I INTRODUCTION............................... 1

Problem Statement .......................... 1
The Piecewise Linear Method................ 4
Illustrative Examples...................... 8
Organization of Dissertation............... 21


An Unpopular Area of Study.................. 23
Phase Plane (2D) Analysis of Friction...... 24
Literature on Piecewise Linear Method...... 25
Extensions to Multi-Dimensional Case
(Approximate Analyses by Describing
Function) ................................ 26
General Results on Periodic Orbits.......... 27
Summary of Literature Review............... 28


System Model with Physical State Variables. 31
Modelling Friction By a Nonlinear,
Two-Input Function .................... 32
The Physical State Model in Piecewise
Linear Regions..... .................... 34
Conversion to Control Canonical Form........ 38
Simplification of DE Solution Via
Coordinate Translation................ 41
Summary of Assumptions and Discussion of
Generality of Model ................... 43



Necessary Conditions for a Simple Limit
Cycle ................................. 49
Illustrative Examples...................... 59
Exact (Necessary and Sufficient) Conditions 68


Stability of Predicted Limit Cycles......... 71
Stability Calculations for Example Problems 91
Symmetry of Limit Cycles for Odd Systems... 100

RESEARCH................................. 116

Results and Conclusions of Current Research 116
Recommendations for Further Research........ 120


LIST OF REFERENCES ................................. 137

BIOGRAPHICAL SKETCH ................................ 139


Figure pa

I-2 TWO-INPUT FRICTION MODEL.............. 3
(B=-l, J=l, LC=1)................... 10
PLOT (B=-l, K=l). .................... 12
ALONG VELOCITY AXIS (B=l).......... 15
(B=1)...... ......................... 16
ALONG POSITION AXIS (B=l).......... 17
(B=.9)..... .......................... 19
ALONG VELOCITY AXIS (B=.9)......... 20

OF F(T1))............................ 65

LOCAL STABILITY NEAR X3............ 96
STATE ............................. 97
STATE (EXPANDED) ....... ............ 98
STATE (EXPANDED) ........ ........... 99

Abstract of Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Doctor of Philosophy




December, 1989

Chairman: Dr. Thomas E. Bullock
Major Department: Electrical Engineering

Exact conditions for the existence of limit cycles due

to friction in a servomechanism have been determined for the

general (n-dimensional) case. These conditions are in the

form of nonlinear algebraic equations, for which a solution

consistent with certain assumptions must be found.

A two-input model for friction is required for accurate

results, including both dry friction and static friction

components. Therefore, standard nonlinear analysis methods

such as absolute stability and describing functions have

limited applicability. A piecewise linear method is used,

which constructs the limit cycle trajectory with pieces from

each region of the state space.

The system model is first converted to a control

canonical form, which simplifies the analysis. Assumptions

made in the formation of the model are specified exactly;

the model represents a rotational servomechanism, but is

shown to be applicable to other systems.


Next, the piecewise linear method is applied to an

assumed limit cycle. The resulting algebraic equations can

be solved for the exact limit cycle period, amplitude, and

trajectory, if any exist. The limit cycles found are of

simple type; that is, complex orbits such as those which

traverse a region more than once would not be found.

Necessity and sufficiency of the existence conditions is


An analysis of the local orbital stability of limit

cycles found by this method is also presented. Eigenvalues

of a given matrix determine the asymptotic stability of the

periodic orbit.

The property of symmetry of limit cycles for odd

systems in general and friction in particular is examined.

Special cases are presented in which symmetry holds, and

counter-examples for the general case of odd systems are

given. The question of symmetry for friction limit cycles

remains open.

Examples are included to illustrate the method. One

example is a 3D system that exhibits two limit cycles in its

phase portrait, each a different type.

Listings of useful analysis computer programs are

given, in addition to suggestions for further research in

this area.



Problem Statement

The objective of this research is to define exact

conditions on the existence and behavior of limit cycles due

to friction. The analysis method to be used is a piecewise

linear approach which, unlike approximate methods such as

describing functions or absolute stability, has the

potential for providing exact (i.e., necessary and

sufficient) conditions for existence.

This problem is of interest in the design of

servomechanisms, which must deal with the nonlinearities in

typical electromechanical devices, such as friction.

Experience has shown that the effect of friction is

generally benign, since it tends to provide increased

damping to the system. However, this cannot be guaranteed

for all systems. Furthermore, it is a standard practice to

precompensate for this effect by deliberately tuning the

system to be underdamped without friction, so that it is not

excessively damped in actual operation with friction. It is

postulated that this technique, when pushed to extremes,

results in large-signal instability, since the damping

effect of friction varies inversely with amplitude. It is,

therefore, a significant research problem to determine the

point where limit cycles can exist.


An approach has been developed during this research

that uses the piecewise linear method, that provides an

exact condition for the existence of a limit cycle in terms

of a nonlinear algebraic system of equations. When a

solution to the system is found (and thus a limit cycle

exists), the method automatically gives exact frequency,

amplitude, and state trajectories for the limit cycle. On

the other hand, exact conditions under which the nonlinear

system of equations is solvable are more difficult because

of the nonlinear nature of the system.

The specific problem under study is shown in block

diagram form in Figure I-1. The plant in the standard

system is represented by a Laplace-transform model of a

motor/load combination, which includes the effect of load

inertia and damping. The model does not include compliance

in the load, although, as described in Chapter III, this may

be included. The inertia (J) and damping (B) represent the

lumped quantities for the load and motor, with any gear

ratio taken into account. The servo controller is

represented as a linear system with motor drive torque as

the output quantity. Note that despite the specific

description of the problem model, it will be shown in

Chapter III (Development of Standard Model) that the model

is actually quite general.

The model used for nonlinear friction includes both

static and sliding coulomb friction components (see Figure

1-2). Therefore, it is a two-input nonlinearity and is



LF= LCsgn[Yn] c

n S
LF= LS sat[LIN/LS ] L

/-- S


thus difficult to handle by standard techniques (such as

describing functions or absolute stability criteria). The

extra complexity caused by the use of the two-input model is

necessary; examples are known where modelling sliding

friction only indicates a limit cycle exists, when, in fact,

the "sticky" effect of the static friction prevents it. On

the other hand, the sliding friction can provide much of the

damping of many servomechanisms which have little or no

viscous friction damping in the load. These comments

indicate that a model including both effects is necessary.

Note that the friction model used is odd, and therefore the

system derivative vector will be an odd function of the


The Piecewise Linear Method

The piecewise linear method is an approach where

systems with piecewise linear nonlinearities can be analyzed

using the tools of linear system theory. Since the

nonlinearity is linear over portions of its domain, the

state space can be divided into regions in which the system

behaves in a purely linear fashion. The linear system

theory then applies, although there will be a different

linear system in each region. Solution trajectories must

then be pasted together from different regions through the

use of boundary conditions (called commutation equations in

relay servo or on-off system analysis such as in Weischedel,


The nonlinear friction to be analyzed in this study is

piecewise linear, even though it requires two inputs. Five

regions of the state space are required (Figure I-3 shows

the case for 2-dimensional state space):

(I) velocity < 0,

(II) velocity > 0,

(III) velocity = 0 and magnitude of accelerating input

torque < LS,

(IV) velocity = 0 and torque > LS, and

(V) velocity = 0 and torque < -Ls,

where LS is the maximum static friction torque (breakaway

torque). The first two regions are always open half-spaces,

while the last three divide the (n-l)-dimensional sub-space

where velocity equals zero into three regions. For example,

for the 3D case (n=3), the last three regions are portions

of the plane y3 = velocity = 0. As can be seen by

examination of the friction model, the nonlinearity either

has a constant output throughout each region (I, II, IV, and

V), or its output is proportional to a linear combination of

states (region III). Even a simple limit cycle must

traverse at least three of these regions, so the linear

systems of each must be considered in this limit cycle


There are several advantages in the use of the

piecewise linear technique. Once a limit cycle is found

(i.e., solutions from four regions that can be pasted

together to form a closed trajectory), the exact trajectory


Y2 = -(K/B) Y1




-LS / K


1 \ I
----- *~~~~~~^






is known. Therefore, the limit cycle amplitude, period, and

characteristics are known exactly, rather than approximately

as in describing functions. Furthermore, the method leads

to a more exact analysis, since the commutation conditions

require the solution trajectories to match exactly at the

region boundaries. Finally, the technique has a geometric

approach that is intuitively appealing and lends itself to

graphical solution.

The method is applicable to the majority of servo

nonlinearities. Backlash, friction, and most other

nonlinearities found in practice are piecewise linear.

The disadvantage of the method is in the difficulty of

solving the commutation conditions for a solution--finding

the trajectories that will "paste" together to form a cycle.

One way to surmount this problem applies the method as a

secondary step in the analysis in order to refine an

approximate solution found by other methods. For example, a

limit cycle period and amplitude are determined

approximately by describing functions; then the piecewise

linear method is used (with the approximate solution as an

initial solution trajectory) in an iterative way to converge

to the true solution. A composite computer analysis program

could integrate these functions for the convenience of the

servo designer. Therefore, the difficulty of solution does

not really limit the applicability of the method.

Illustrative Examples

Before going into the detailed development of the

method and other results, it is useful to present simple

example systems exhibiting friction limit cycles. The limit

cycles were analyzed by the methods to be presented later,

so we can return to these examples at the appropriate time

to illustrate the method; here only the results are shown.

The second example system is particularly interesting

since it is only 3D (three state variables describe the

system), yet it exhibits two distinct limit cycles of

differing types (when the appropriate parameter values are

selected). The fact that such a simple system can have such

complex behavior indicates the richness of the study of

piecewise linear systems.

Example I-1: Two-Dimensional (2D) System with Sliding

Friction Limit Cycle

The question of friction limit cycles for the case n=2

can be completely solved; in fact, the piecewise linear

method is the same as standard phase plane methods for this

case. For n=2, the controller can be a constant gain only

(since there are two states in the plant); the state

variable model of the system is then

(1.1) dyl/dt = y2

dy2/dt = -(K/J) yl (B/J) y2 LF/J

where LF is the friction torque from the model in Figure

1-2, and any feedback gain on y2 is lumped in with B. It

can be shown (see Thaler and Pastel, 1962, pp. 96-104) that

for B>0 and for K>0, the system is globally asymptotically

stable (i.e., no limit cycles). Thus a stable system

remains stable (in fact it is more stable) with friction (2D

case only).

The case B=0 is also stable; the linear system has

imaginary eigenvalues, the friction gives enough damping to

give asymptotic stability (trajectories spiral to an

equilibrium point that is not necessarily the origin).

Interestingly, a limit cycle exists for every case with

negative damping (B<0). If the position feedback gain K is

sufficiently large to give the linear system a complex pair

of poles, an unstable, unique limit cycle exists, whose

amplitude depends on K, B and the inertia J and sliding

friction level LC. This "sliding" type of friction limit

cycle is one in which the motor or other object being moved

never "sticks"; the reversing torque is sufficient when it

comes to rest to immediately breakaway in the opposite

direction. Therefore, the orbit consists of two parts, each

representing motion, with two switching per cycle. The

initial conditions can be expressed in terms of the system

parameters as

(1.2) X10 = (Lc/K) [(L+1)/(X-l)], x20 = 0


(1.3) p = exp[-Br/(2JP)], 3 = [4K/J B2/J2]k

Figure I-4 shows some values of amplitude (which equals x10)

as a function of gain K for the case J = LC = 1, B = -1; it





S-- Q


oo II

u- E- II




__________________g 0

O a

appears that a limit cycle of any amplitude can be found by

adjusting gain properly.

The unstable limit cycle bounds the region of

asymptotic stability: trajectories inside spiral in to

equilibrium near the origin, while those outside go to

infinity (Figure I-5 shows phase plane plot). This example

shows why it is dangerous to design underdamped systems, and

depend on the sliding friction to stabilize the system; the

friction damping is amplitude-dependent, and at some

critical amplitude, is insufficient for stability.

Example 1-2: Three-Dimensional (3D) System with Both

Sliding and Sticking Limit Cycles

The second example system to be considered is defined

by the differential equations

(1.4) dy/dt = A + b LF


(1.5) A = 0 1 0

0 0 1

-K1 -K2 -B

with b' = [ 0 0 -1 ], LF = friction torque, and where the

state variables are y2 = position, y3 = velocity, and

yl = compensator integrator. Figure I-6 shows a block

diagram of this system.

Case I: One Stable Pole and One Imaginary Pole Pair (Only

Sticking Limit Cycle Exists)

For specific results, let us first set B = K1 = K2 = 1,

and set the sliding friction torque to 1.0 and sticky

o o 0 0 0
S0 C

> O OU u






c I



oi a
3 &<



friction (breakaway torque) to 1.2. Simulation of this case

demonstrates a sticking limit cycle with an initial

condition of y0' = [ 1.0 0.2 0.0 ]. This type of friction

limit cycle repeats a four-part cycle of sticking (until

torque integrates up to breakaway), sliding to a new

position, sticking again, and sliding back to the start

position. A visual display of this behavior can be seen in

Figures I-7 through 1-9, which show projected views of the

closely similar limit cycle found using the Case II


The solution of the piecewise-linear equations results

in a valid limit cycle solution where y = [ 1 0.2 0 ],

T1 (the sliding period) = r, and the sticking period

T2 = 10, yielding an overall limit cycle period of approx.

26.3 seconds, matching the simulation results.

Case II: One Stable Pole and One Unstable Pole Pair (Both

Sticking and Sliding Limit Cycles Exist)

If the same example system is used with B = 0.9, and

K = K2 = 1 still, the closed-loop eigenvalues (of the

linear portion of the system) move into the right-half

plane, to .026 j 1.024, while the other pole is at -.9524.

As might be expected, there is still a sticking limit cycle

solution close to that of Case I (Figures I-7 through I-9).

Numerical methods give a solution at approximately

O = [0.98 0.22 0 ], and approximately the same period

(confirmed by simulation).

1- + +


I - I








i< H









______ _________ __________ _______ ----





> c OU-OU E


C E-

1 0







CC 3

___ __ __ _ _ _ _ _ _ _4_____ I__ _ _

t t + +



o o





The interesting point about this case, however, is the

similarity to Example I-1, the 2D case; that is, an unstable

complex pole pair exists. Since the nonlinear damping

effect of the sliding friction is largest at small

amplitudes, we might expect trajectories to converge for

small amplitudes and diverge at large amplitudes where the

unstable linear poles overcome the nonlinear damping (as in

the 2D case). Therefore, this case would have both a

sticking and a pure sliding friction limit cycle!

Simulation having shown a limit cycle with a

half-period T1 approximately equal to 3.1 seconds, numerical

solution of the equations for the symmetric, sliding case

was performed, yielding a solution at T1 = 3.1455 seconds,

and y0 = [-0.283 12.4 0.]. Figures 1-10 and I-11 show

some views of the two 3D limit cycles for this case.

It is expected from physical intuition and simulation

results that the sticking orbit would be stable, while the

sliding orbit would be unstable. This is in fact the case,

with the sliding orbit exhibiting saddle point behavior--one

stable and one unstable mode (real eigenvalue / eigenvector

pair), while the sticking orbit has stable node behavior--

two real stable modes. The analysis of the orbital

stability will be examined in Chapter V, where Example 1-2,

case 2 will be used to demonstrate the method used for limit

cycle orbital stability calculations.




S -1

-- O

- o
-- 1 I -- ---- -- \ I -- -- I -- I I I ---- -- -- I \ ---- -

> OU)HE-q >





- 4


-- 0

0 0 0 0
-4 -4 ( )








Q &>

"Z=E-HuC < E-,OW OE-4

Organization of Dissertation

The main body of this dissertation is organized as

follows. After a review of existing results on this problem

and related areas, the standard model introduced above is

developed more fully. The objective in the model

development (in Chapter III) is to provide more convenient

forms for analysis (a control canonical form and a normal

form), while maintaining the generality and traceability to

actual servomechanisms of the original model.

Once the model is in the proper form, Chapter IV

develops the nonlinear system of equations that represent

the necessary and sufficient conditions for the existence of

a friction limit cycle. Several special cases are examined

that are somewhat simplified from the general problem and

easier to solve. For example, if it can be assumed that the

limit cycle is symmetric about the origin (that is,

half-wave symmetric), the equation system can be simplified


The benefits of symmetry in the solution of the limit

cycle conditions led to an examination of the conditions

under which this symmetry existed. Chapter V examines this

problem and describes special cases for which symmetry

exists. Unfortunately, it is not known if friction limit

cycles are symmetric in general, although a search for

counterexamples was unsuccessful.

Using the results of the previous chapters, a complete

stability analysis of friction limit cycles is also included

in Chapter V. Using the eigenvalues of the local stability

matrix developed in this analysis, the stability of limit

cycles previously found is determined. This allows the

global phase portrait of the system to be pieced together

from a definition of all the equilibrium points, limit

cycles, and their stability characteristics.

The results presented here are considered only a start

in an essentially new area: use of piecewise linear

techniques for limit cycle analysis. Many possibilities for

further research exist, both in the completion of the

solution for friction, and in application of the method to

other piecewise linear nonlinearities. As pointed out

above, the majority of the plant nonlinearities encountered

in servo design are piecewise linear, including gear

backlash, saturation, quantization, and dead zone. The last

chapters provide a summary of results and suggestions for

further research in these areas.


An Unpopular Area of Study

Surprisingly little work has been done in this area,

for various reasons. First, the nonlinearity is not "nice,"

i.e., continuous or smooth, which prevents the application

of many mathematical tools. Second, when modelled with two

inputs it really becomes a system with multiple

nonlinearities, for which standard approaches cannot be

used, and extremely few results are available.

Finally, and, I believe, most importantly, it is

perceived as generally benign in effect. It can be shown

fairly easily (see Thaler and Pastel, 1962, for example) in

the 2D case (i.e., two state variables), that if the linear

portion of the system is already stable, the system with

pure sliding friction added is also stable. .It seems

probable that this result could be extended to arbitrary-

dimensioned systems, based on the following energy argument:

the linear system is stable (hence contracting), and because

the (sliding) friction opposes motion, it damps the system

further, so the map is more contracting than without it.

Therefore, friction is generally held to be benign in effect

(it increases the stability of a system by increasing the


For these and possibly other reasons, the references

available in this area were few, and no analysis of the

general multi-dimensional friction limit cycle problem was


In spite of these comments, the study of friction limit

cycles has practical value for control system design.

Although the effect of sliding friction is to increase

damping, it has become fairly standard to precompensate for

this effect in a design. The system is deliberately

designed to be underdamped, so that it is not too well-

damped when friction is added. It is conjectured that this

technique, when pushed to extremes, results in large-signal

instability, since the damping effect varies inversely with


In addition, the stiction component can cause

difficulties also. An example was given in Chapter I of a

stable (although not asymptotically stable) third-order

linear system which has a limit cycle when nonlinear

friction (including stiction effects) is added. A slight

variation in parameters results in a stable linear system

which destabilizes with the addition of nonlinear friction!

Phase Plane (2D) Analysis of Friction

As mentioned above, Thaler and Pastel (1962), in their

classic text on nonlinear systems, completely solved the

second-order case for friction, including both sliding and

static components. They also give an exact criterion for

the existence of limit cycles when an input ramp is present.

They showed that the system had no limit cycles (for a

system with positive damping) for the zero-input case.

Earlier work by the same authors, (Pastel and Thaler, 1960)

actually shows a stability boundary for backlash (as a

function of system damping), then demonstrates the

stabilizing effect of coulomb friction.

Analyses of the 2D friction problem were presented in

several other references from the period when phase plane

analysis was an active area of research.

Literature on Piecewise Linear Method

The approach used in this dissertation to give exact

solutions for limit cycles (if they exist) has been known

and used for many years. In the relay servomechanism

literature, it is known as the piecewise linear method.

The method was used primarily for relay systems, though

it is generally applicable to any piecewise linear

nonlinearity. Note that the piecewise linear method is

essentially the same as phase plane analysis for 2D systems,

so there is overlap between the phase plane references, such

as those cited previously, and those for the piecewise

linear method. Unfortunately, the graphical methods of

phase plane analysis cannot be easily extended to higher

dimensions (see comment below on the work of Ku).

The piecewise linear method was used at least as early

as 1963 (by Kovatch); additional studies can be found in

references such as Kovatch (1964) for two nonlinearities,

O'Donnell (1964) for time-optimal switching, Marstrander

(1969) for backlash, Negoescu and Sebastion (1971), Langill

(1965), and Urabe (1967). Weischedel (1973) applied the

method to on-off control systems. Ku (1958) made

significant and early contributions, including attempts to

extend the phase plane method to higher dimensions by

various graphical projections.

The Russian controls literature calls this the Andronov

point transformation method; a reference was found in a

translation of material developed in 1956 by E.P. Popov

(translated 1962). Although significant development of this

area has been performed, especially in the Russian journal

Automatika, an application to the friction problem has not

been found.

Extensions to Multi-Dimensional Case
(Approximate Analyses by Describing Function)

Analysis using approximate techniques such as the

describing function go back at least as far as the piecewise

linear method, or even farther in the Russian journals

(Popov, 1956). Four papers are listed by Gibson (1963),

that he refers to as independent developments of the

describing function method; all are from the late 1940s!

However, the two input difficulty that limits the

applicability of describing functions to friction is usually

dealt with by ignoring the static friction. No references

have been found that discuss the general case for an exact

friction model. No other references were found that

attempted the multi-dimensional friction problem with other


General Results on Periodic Orbits

In the mathematics literature on the theory of periodic

solutions of ordinary differential equations, the Poincare

map (Hirsch and Smale, 1974) is frequently used; this map

represents the change in system state after one cycle or

period. By searching for fixed points of the map (which

represent a return to the original state after one cycle), a

periodic solution can be found. This theory was used in

generating the stability results to be stated in Chapter V,

where the orbital stability of friction limit cycles is

examined. The Chapter V results are, therefore, original

only in their application to friction limit cycles, not in

the theory of orbital stability.

There are many results in this literature (although the

Poincare theory itself is restricted to planar systems).

Pliss (1966) stated some results on the existence and

uniqueness of these fixed points, although they could not

apply to friction as they were restricted to continuous

functions, or non-autonomous systems. Hale (1969) discussed

fixed-point theorems for 2D systems and n-dimensional

systems on a torus; in addition he stated the orbital

stability result (that all eigenvalues of the map except one

must be less than one in magnitude). Hayashi (1964) refers

to this approach as the topological approach to periodic

solutions, since the theorist is concerned with the

existence and topological properties of integral curves in

state-space. Finally, Hirsch and Smale (1974) provide an

excellent discussion on an introductory level of the theory

of flow maps; this reference was invaluable in the stability

analysis of Chapter V.

My survey of this literature was thorough, but not

exhaustive; no applications to this nonlinearity were


Summary of Literature Review

In the 20 years since the heyday of phase plane and

piecewise linear methods, it appears that almost no new

results have been discovered in this area. However, the

study of friction limit cycles has significant practical

value for designers of servomechanisms, since a servo

designed using standard techniques may exhibit instabilities

in the presence of friction.

The literature survey was extensive in the area of

controls, and moderate in the mathematics field. In the

former, IEEE indices were searched, including 25 years of

the Transactions on Automatic Control. Twenty-five years of

the Russian journal Automation and Remote Control were

examined. The journal Automatica was also reviewed.

Although no journals in the field of pure mathematics were

checked, many texts in the areas of stability, limit cycles,

and dynamical systems were examined.

The main contribution of this research is the

application of known (although relatively obscure in

controls literature) techniques such as the piecewise linear

method and Poincare flow maps to the problem of nonlinear

friction. No reference was found that considers the multi-

dimensional friction limit cycle problem using an exact

method such as the piecewise linear technique. All the

references found either analyze the 2D case for friction, or

derive approximate results by using a simpler friction

model. Practical experience with servos indicates that it

is unacceptable to ignore either friction component.

The development presented here, therefore, appears to

be the first to attack the multi-dimensional friction limit

cycle problem with sufficient accuracy for practical

applications. In addition, the results on limit cycle

symmetry for piecewise linear systems, and the derivation of

limit cycle orbital stability by the same method, are

original, to the author's knowledge.


This chapter is concerned with setting up the model of

the system to be studied, and converting it to a convenient

form for analysis. The original system model is in the form

of physical state variables, where the states represent

actual physical quantities in the servomechanism. The model

form to be used for most results is a control canonical

representation, which is developed from the physical states

by a linear transformation. This form has advantages in

clarifying subsequent developments of the theory, including

the role of system zeroes and characterizing limit cycle

symmetry. This chapter also demonstrates how the piecewise

linear friction torque input can be removed by a

translation, which leads to a simplification of the

existence formulas of Chapter IV.

In addition to the derivation of the state forms, this

chapter provides a description of the various assumptions

made in the analysis. In other words, the generality of the

analysis and the solutions are defined by examining the

assumptions made. For example, although the original

motor/load model included no compliance or resonance

effects, this restriction can be removed, allowing the

results to be more generally applicable.

System Model with Physical State Variables

The two primary physical variables are angular velocity

(required as input to the sliding mode friction) and its

integral, angular position. The source of this model is a

rotating servomechanism. The theory is not restricted to

rotational friction problems, however, since angular

variables (torque, angular velocity, etc.) could be replaced

by linear ones (force, velocity, etc.) without affecting the

nature of the problem (this point is discussed further at

the end of this chapter). The model was illustrated in

block diagram form in Figure I-1; the controller transfer

function G(s) is replaced by state equations in the model


The remaining n-2 state variables, where n is the order

of the system, are controller states. That is, they form

the model of the servo loop compensation, amplifiers, motor

armature effects, and other dynamic portions of the system.

Thus we may set up a model of the following form:

(3.1) = Ac y + a y + a c2yn
c c -cl n-l -c2 n

Yn-1 = Yn

J Yn = -Kyn-1 Byn + b'yc + f(y)

where the yi variables are the physical states and make up

the vector y, the dot (') indicates the time derivative

d/dt, prime (') indicates transpose, f(.) is the friction

torque, J is the lumped system inertia, and B is the

damping. The n-2 dimensional linear subsystem in the n-2

states in vector yc is formed by the system matrix Ac
c c

(n-2 x n-2), the input vectors a cand a (n-2 x 1), and

the output vector b (n-2 x 1) that feeds the output of the

controller into the torque equation.

The only assumptions made in the construction of this

model are that yn-1 and yn are angular position and

velocity, respectively. These assumptions completely define

the first equation, and require that the nonlinear friction

torque term f/J appears in the y differential equation

(representing angular acceleration due to the sum of torques

applied). The rest of the system is linear; the friction

term does not appear in any other derivative.

Note that there is an assumption made here that

restricts the generality of the model: that the friction

appears in only one equation. Although examples can be

constructed that do not obey this restriction (such as those

involving differentiators, full state feedback, torque

measurement, or a foundation model), this model is felt to

be quite general, and covers essentially all servo problems

where friction limit cycle information is of.interest. This

point will be discussed in more detail later in this


Modelling Friction By a Nonlinear, Two-Input Function

The nonlinearity investigated in this analysis is a

model for mechanical friction that depends on both the

relative velocity of the surfaces and the driving torque of

the moving element. The level of friction opposing the

motion depends on velocity (or more precisely, the direction

of the velocity) when the element is moving. In the case

where the element is stationary, on the other hand, the

static friction torque holding the system motionless must be

set equal to the input driving torque that is attempting to

move the element. Each of these cases is based on the

models of Coulomb friction described in any basic statics

textbook. Figure I-2 in Chapter I illustrates the friction


Although an analysis can be performed with only one

element modelled, thereby using a single-input nonlinearity

and simplifying the analysis, this would lead to results

that would not apply in any physical system, and could be

misleading or erroneous. Cases can be constructed where the

modelling of sliding friction alone predicts a limit cycle,

due to its highly nonlinear, discontinuous characteristic.

In a practical system, however, no limit cycle would exist.

For example, the two-dimensional friction limit cycle

is completely solved in Thaler and Pastel (1962), as

described in Chapter I. Neglecting the static friction

component in this analysis leads to the conclusion that a

limit cycle of very small amplitude can be obtained.

However, in reality, such a small initial condition will

stick at the initial point due to the static friction


On the other hand, the sliding component of friction

cannot be ignored either. As a qualitative example of a

case demonstrating this, consider a system with no damping,

which would limit cycle if static friction alone were

present. The system starts at an initial position with

nonzero error, integrates up in the controller until it

breaks free from the static friction, and overshoots the

desired position to start the cycle all over again. The

inclusion of sliding friction in the model, however, might

provide enough damping to the system to cause it to spiral

in to the equilibrium point at the origin. Another

possibility indicated by the behavior of typical friction

limit cycles found is that a system that is actually stable

for small signals (bounded region of stability) might appear

unstable if analyzed without sliding friction.

These cases indicate that accurately predicting limit

cycles for physical servomechanisms requires a complete

model of friction, as is used in the analysis contained in

this paper.

The Physical State Model in Piecewise Linear Regions

As described in the introduction, application of the

piecewise linear method to the system model results in five

regions of state space, with a linear system for each

region. The form of the (now linear) model in each region

is described, so that the trajectories in each region can be


Region I: Velocity < 0 (v <0). This n-dimensional region is

described by the same set of differential equations as in

(3.1), except that the nonlinear friction term f(y) is

replaced by a constant input Lc. This follows from the fact

that the friction in sliding mode is constant (independent

of velocity) as long as the system is in motion in one

direction. From basic mechanics, this constant friction is

a function of the normal force, with the proportionality

constant being the coefficient of sliding friction. The

torque Le is positive, since physical friction always acts

to oppose the direction of motion; velocity is negative, so

the acceleration torque due to friction is positive.

Therefore, the model in region I is a linear system driven

by a constant input.

Region II: Velocity > 0 (y >0). This region model is

exactly the same as in region I, except that f(y) is

replaced by -Lc

Region III: Velocity = 0 (y =0) and Acceleration Torque 5

Breakaway Torque. In this case, the system is in a static

condition, and the model for the static component of

friction applies (also known as stiction). The friction

opposes the applied acceleration torque with an equal torque

so that the net acceleration torque is zero,-up to a maximum

amount of static friction. The maximum amount of friction

(the breakaway torque) is determined by the coefficient of

static friction and the normal force. Once the acceleration

torque becomes greater than the maximum stiction, the net

torque becomes nonzero, and the system begins to move

(enters region I or II).

Note that the acceleration torque used to define this

region, and which is matched by an opposing friction force,

consists of every term in the torque equation (the

derivative of y ) except the friction f. Therefore, the

linear terms in the torque equation define the boundaries of

this region. The region is contained in the

(n-l)-dimensional hyperplane y n=0.

Since the net torque is zero and the velocity is zero,

the first two differential equations in (3.1) drop out (the

state yn-i remains constant and y remains zero over

trajectories in this region), leaving the n-2 dimensional

controller system as the only equations required to define

the trajectory in region III. Of course, input vector ac2

drops out (since yn=O), while acl and b are used to model

the effect of the constant angular position on the

controller and the controller output torque (to determine

the breakaway condition).

Region IV: Velocity = 0 (y =0) and Acceleration Torque >

(+)Breakaway Torque. This region is also contained in the

(n-l)-dimensional hyperplane Yn=0, but the driving

acceleration is sufficient at any point to overcome the

maximum static friction opposing the impending motion. Thus

a system that comes to rest with sufficient torque

(trajectory enters region IV rather than III) immediately

breaks into motion again, entering region II (since

acceleration torque is positive, velocity increases from

zero and trajectory must enter region II).

For this reason, a trajectory can be said to cross this

region (in the sense of crossing a plane, for example), but

is only in the region for an instant. Since the integral of

a finite quantity over a set of measure zero is zero, the

applied torque and the differential equations are irrelevant

to the system behavior. The trajectory leaves the region at

the exact point it entered. It is therefore unnecessary to

examine the form of the system equations in this region.

Region V: Velocity = 0 (y =0) and Acceleration Torque <

(-) Breakaway Torque. The region V model is identical to

the region IV model, except for the breakaway conditions.

The five regions cover the state space, including the

equilibrium point (actually, line segment) at the origin.

The trajectories in each of the individual regions can be

described by the solution to the linear system of

differential equations:

(3.2) y(t) = F(t,t0) y(t0) + G(t,t0) L

F(t,t0) = exp[A(t-t0)]

where L is the (constant) input in that region, A is the

system matrix for that region, and the F and G matrices vary

from region to region. The initial conditions y(t0 ) are

determined by the commutation conditions at the boundary of

the region at which the trajectory enters.

An algebraic trick, to be presented below, can be used

with this form to remove the input term (GL) in equation

(3.2). The system then behaves as an autonomous system,

with considerable simplification of the defining equations.

Conversion to Control Canonical Form

The most convenient form for analysis is the control

canonical form, because the role of system poles and zeroes

in the existence conditions in Chapter IV will be clarified

through the use of this representation. In addition, a

result in Chapter V on symmetry is proved using this form.

Although this section constitutes a digression from the main

line of development, it is justified by the usefulness of

the control form.

The conversion to control form is accomplished by

representing the controller (transfer function G(s) in

Figure I-1) in control form and combining this

(n-2)-dimensional system with the equations for position and

velocity. When completed, the system model will have the

following form:

(3.3) x1 = x2

X2 = x3

Xn-1 = xn

xn = -a0X1 alx2 .. axn + f()

The coefficients in the differential equation for x are the
coefficients of the closed loop polynomial of the linear


The procedure used to convert the system to control

canonical form is as follows. The physical state model is

set up so that the controller (the (n-2)-dimensional

subsystem represented by state variables yc) is in control

canonical form already. This can be done for any controller

representation that is controllable (if the original

representation for the controller was a transfer function, a

controllable representation can always be found).

Next, the model is transformed so the velocity feedback

into the controller subsystem is zero (i.e., vector a in
(3.1) is zero). This can be done without loss of generality

in the model, as shown by the following argument: Suppose

ac2 is not zero in the original system. We can apply a

similarity transformation

(3.4) y = Ta (1)

where y(1) is the new state variable vector, and

T = 1 0 0 . 0 a 0

0 1 0 Oa 0

O0 0 O l a 2 0
0 0 0 1 0

0 0 1

Note this is equivalent to replacing all (n-2) controller

states yc by (Yc(1) + ac2 Yn-l)- A straight forward

calculation (plug (3.4) into (3.1)) shows that the resulting

system has no feedback of velocity into the controller.

Since a similarity transformation (or the equivalent

change of variables) does not affect the basic system

behavior, but merely its representational model, this

alteration of the model is valid. Physically, this

operation has removed a redundant parameter in the original

model, since the velocity feedback in the original model can

be achieved by adjusting the damping parameter, B, or

adjusting the position feedback gain, K1. Since the models

are equivalent, the superscript on the y state variables can

then be dropped in subsequent equations.

With the controller in canonical form, and the velocity

feedback to the controller eliminated, the following

intermediate form is obtained:

(3.5) yc = Ac + a yn-

Yn-1 = Yn
J n = -K1n-1 By + b'c + f()

In order to complete the conversion of the model to

control form, two similarity transformations are then

applied to this intermediate model in succession. As

previously explained, similarity tranformations do not

change the basic behavior, hence are valid alterations of

the model equation. The transformations are:
(3.6a) A, = T1 A T1
(3.6b) A2 = T2-1 A T2

where the nxn matrix A is obtained from the model in (3.5),


(3.7a) T1 = 1 0 0 0

0 1 0 0

0 0 1 0 0

0 1 n-3 1 0
0 .0 0 1

(3.7b) T2 = 1 0 0 0

0 1 0 0

0 0 1 0

0 60 1 n-3 1

and where the 's are the coefficients of the

characteristic polynomial of the controller alone (i.e., the

characteristic polynomial of matrix Ac in equation 3.5).

The transformed system is then in control canonical form.

Note that the friction input is unchanged by these

transformations, and is still applied in the yn differential

equation. The composite transformation is required for

later derivations, as it defines the relationship between

the physical states and the control form states:

(3.8) y = T1T2x


T = TT = 1 0. 0 0

0 1 0 0

0 1 0 0

P0 16 n-3 1 0
0 60 91 n-3 1

Simplification of DE Solution Via Coordinate Translation

Once the system is in control form, the following

translation of the state coordinates can be used to convert

it to an autonomous system (within each region). Since the

friction is a constant in region I or II (L or -L ), and
since none of the differential equations except that for x
involves xl, the xl coordinate can be translated by

(3.9) xla = x1 f/a0

where the subscript "a" indicates autonomous coordinates,

and f is the (constant) friction level in that region (I or

II). Substitution of this translation shows that the x
equation still gives the same result, the system matrix is

unchanged, and the input has disappeared. The autonomous

system model is then

(3.10) xla = x2

k2 = 3

n-1 = n

x = -ax x nx
n 0 la 1 2 n n
Therefore, by performing a translation on the initial

condition and on the final state of the trajectory in the

region, the trajectory can now be described by the simpler


(3.11) x (t) = F(t,t0) a (t0)

F(t,t0) = exp[A(t-t0)]

instead of equation (3.2), where the matrix A is in control

form (i.e., the matrix A2 produced by transformation 3.6).

Of course, the control form representation is not

required to perform this translation. The equivalent

translation in terms of the physical state variables can be

defined using the transformation matrix as

(3.12) ly = T2 1

Ya = Y ly

where 1' = [ f/a0 0 . 0 ] is the translation of the

control form states (equation (3.9)) and Ya is the

autonomous physical state variables. The resulting

differential equation for y would be of the form

(3.13) ya = A Ya

(with A unchanged). A slightly different form of the

existence conditions of Chapter IV would then result. The

control form has certain advantages, however, as will be


Summary of Assumptions and Discussion of Generality of Model

The following assumptions and conditions on the system

are made during the analysis of friction limit cycles in

this dissertation:

(Al) The system can be modelled by a model of the form

in Figure I-1, and equations (3.1).

(A2) The velocity feedback into the controller is zero

(without loss of generality, as shown by argument given

previously in this chapter).

(A3) The controller portion has a controllable

representation (this is assured if the controller is

representable by a strictly proper transfer function; if

transfer function is proper the direct feed-through term can

be divided out and lumped with the direct position feedback,

so this case is also allowed).

(A4) Friction appears only in the differential

equation for velocity (the acceleration torque equation).

(A5) The closed-loop, linear portion of the system

should not have a pole at the origin.

(A6) The limit cycles to be analyzed are all simple,

that is, they traverse regions I and II once before

returning to the initial point, and therefore have only four

switching per cycle at most.

Note that many of the assumptions were made for

convenience in the subsequent development, so that a

specific model could be used for concreteness. Most can be

removed by a transformation of the original problem as

discussed below, or by minor modifications of the theory,

and so involve no loss of generality.

In fact, the only restriction to the system model that

cannot be removed is that the limit cycle is simple

(assumption A6). The analysis could be extended, as

discussed below, to remove even this restriction. Thus the

theory is quite general.

Several of these assumptions require further

discussion. Assumption (A4) is required to perform the

analysis using the piecewise linear method in Chapter IV,

since the analysis uses the control canonical form heavily.

If the friction input appears in more than one differential

equation, the trick used there to eliminate the friction by

a coordinate translation is invalid.

Examples could be generated that violate this

assumption. If acceleration torque was directly measured,

the friction torque could be inferred and used as an input

to the controller. The same comment holds for systems with

differentiators, since the derivative of the velocity

involves the friction. Friction also causes a reaction

torque on the mounting of the motor or other physical

device, so if a foundation or mount model is included, the

assumption would be violated.

This condition does not restrict the generality of the

method, however. Systems which violate this assumption can

still be transformed to control form, the method applied,

and existence conditions checked. A minor modification of

the Chapter IV conditions (involving the P vector only)

would be required.

Assumption (A5) is not necessary, except for steps

involving taking the inverse of certain matrix quantities in

Chapter IV. However, the assumption results in no loss of

generality, for the following reason: a system with a zero

eigenvalue has a free integrator so that, with the proper

selection of state variables, one state does not feed back

into the system. This can be seen by examining the control

form of the system model; when an eigenvalue is zero, the

coefficient a0 of the characteristic polynomial is zero, and

the first column of the system matrix is zero (see equation

3.3). Therefore the magnitude of the state xl has no effect

on the behavior of the system. Any limit cycle existing in

the system will thus also exist in the (n-l)-dimensional

system formed by deleting the state xl (although not

necessarily vice versa). Therefore, any system that fails

this assumption can be analyzed by deleting free integrator

states, then looking for limit cycles; any limit cycles

existing in the original system would be found by this

method. Note, however, that solutions may be introduced

that do not exist in the original. For example, an

asymmetric limit cycle in the reduced system might not exist

in the original since the free integrator output might not

be periodic (see Chapter V).

Assumption (A6) is also necessary to allow the

derivation of Chapter IV to proceed. The piecewise linear

method requires the overall shape of the limit cycle to be

known, including which regions are traversed and the order

of traversal. This could be considered a limitation of the

method, since the analysis cannot be performed without

knowing the shape of limit cycles that are to be found by

the analysis!

Note, however, that an even more restrictive assumption

is made in the classical describing function analysis, a

tool that has proven quite useful for decades. That

analysis assumes a sinusoidal input to the nonlinearity,

therefore implicitly assuming a half-wave symmetric type of

solution with two switching per cycle. It is therefore

felt that the assumption is reasonable. Future research

could remove this assumption by listing all possible

solution shapes (and order and number of traversals of

regions during each cycle), and applying the method to each

possibility. The difficulty lies in the large number of

potential solutions to be checked.

If it could be proved that all limit cycles are simple,

the difficulty would disappear. Without such a proof, the

results of this thesis apply only to simple limit cycles.

For example, the necessary conditions for the existence of a

limit cycle found in the next chapter are necessary only for

simple limit cycles. Therefore, the lack of a solution for

these conditions for a particular system only guarantees

that a simple limit cycle cannot exist for that system.

More complicated behaviors are not ruled out.

Finally, some comments on the generality of the system

model assumed are in order (assumption (Al)). The model was

inspired by the case of a rotational servomechanism, such as

a DC motor. Although this is an important application, the

theory developed here can be applied to many other systems,

as long as they can be put in the form specified. For

example, a linear drive system has the same characteristics

of inertia, damping, friction, and acceleration force

(instead of torque), so the theory can be directly applied.

Systems with gear-trains and loads can be modelled by

lumping the load and gear inertia into the motor, as long as

the system is rigid and closely coupled (gear backlash can

be ignored). Although the model is for a continuous-time

system, the use of discrete-time (i.e., digital) controllers

can be accomodated by using their continous-time


Many systems of practical interest have mechanical

flexure, modelled by resonances. If they can be lumped into

the model of the controller, the theory can still be applied

as is. However, the case of two masses, separated by a

spring, where both masses experience nonlinear friction,

does not fit into the present model. Additional piecewise

linear regions must be specified, with a great increase in

the complexity of the analysis, unless the limit cycle is

known to be low in frequency so the system could be

considered to behave as a rigid body.

In summary, certain assumptions are made to allow a

specific, concrete model to be used in the development.

Examination of the assumptions shows that most are for

convenience only, and the theory applies to a quite general

class of problems.

Given these restrictions and assumptions on the problem

to be studied, and the development of a convenient

representation, the piecewise linear method can be applied

to the problem of friction limit cycles. The next chapter

derives exact conditions for the existence of these limit

cycles, and, once a solution is found, provides exact

information on the limit cycle trajectory.


This chapter applies the piecewise linear method to the

problem of friction limit cycles, using the model developed

in the previous chapter. The method provides a set of

nonlinear algebraic equations for which a solution must

exist in order to have a limit cycle. These represent

necessary conditions for limit cycle existence, since every

(simple) friction limit cycle for this system model meets

these conditions. The issues of sufficient conditions and

minimal sets of necessary and sufficient conditions are

discussed. Finally, some examples applying these conditions

are presented.

Necessary Conditions for a Simple Limit Cycle

In order to simplify the analysis, the assumption is

made that the limit cycle which the system exhibits is

simple, in the sense that it traverses each of the four

state-space regions exactly once during each limit cycle

period. This is assumption (A6) discussed in Chapter III.

This assumption makes sense from a physical viewpoint, since

it is hard to imagine a system where the velocity changes

sign more than twice before returning to the initial state.

It is clear that this assumption holds in the 2D case,

since the trajectory must go to the left in region I and to

the right in region II (see Figure I-3 and I-5) and the


trajectory cannot cross over itself. However, no claim can

be made that the following analysis is completely general,

due to this assumption. Therefore, although this is a

restriction of the generality of the results, the solution

of the equations requires the limit cycle meet the following


Simple Limit Cycle: A periodic orbit which traverses each

piecewise linear region once at most (at most twice for

region III).

As a consequence of this assumption, there will be

exactly four distinct switching instants during the period.

The system will be in a particular region, and, therefore,

its behavior will be governed by that region's system model

equations, for T. seconds, for i = 1, 2, 3, 4. As will be

seen below, some of these switching periods may be zero.

Thus the total limit cycle period T is the sum of the four

(4.1) T = T1 + T2 + T3 + T4

Note that the subscript on the periods indicate order, and

not the region being traversed during period i.

Since the limit cycle is periodic, any initial point

may be chosen to start the analysis. Define x to be the

system state at t=0. Assuming a limit cycle exists, this

initial state can be chosen to be the point in region III

(or V, see Figure 1-3, p. 6) where motion is impending. In

other words, initial velocity is zero, and the initial

driving torque is negative, and sufficient to overcome

static friction and begin motion. These requirements on the

initial point can be expressed as (recalling that y

represents the physical state variables, and x the control

canonical states):

(4.2a) Yn(0) = [ 0 60 n-3 1 ] x

= 1' 0X = 0

(4.2b) yn(0) = H' A x" < -LS (or = -LS )

where the notation A' is used to represent the nth row of

the x to y transformation in (3.8) (used to calculate the

velocity yn from the state x), the parameter LS is the

maximum static friction before breakaway, and the matrix A

is the control canonical form of the state matrix, in region

I. Note that the A matrix from region I (velocity < 0) is

used since the rate of change of yn is identically zero in

region III, by design.

Having characterized the requirements on the initial

point, the limit cycle trajectory can be followed in region

I (which the system transits after breaking free) using the

region I model. The autonomous state variables can be

formed as discussed in Chapter III, by making the


(4.3) xla = x (-L/an)


x = x -

1' = [ (-Lc/an) 0 0 .

0 ]

where a is the value in the nth row and first column of A,

and LC > 0 is the sliding coulomb friction torque. This

translation is selected so that the friction torque is

positive, and opposes the motion in region I in the negative

direction. Note that the parameter anl cannot be zero,

since this would give a zero column in A, which signals the

presence of a free integrator.

The trajectory in region I is then
At At
(4.4) x (t) = e x = e ( 1)
-a aO (
so that the trajectory exits region I at the point where y

(velocity) again goes to zero

(4.5) x(T1) = xa(TI) + 1
= e (x 1) +1

with the requirement that n (T ) = 0:

(4.6) n(T1) = l' X(T1) = 0

At this point, the limit cycle will display two types

of behavior, depending on the driving torque at t=T If

the trajectory transits from region I into region V, this

torque is already sufficient to break loose the system and

drive back in the other direction, then T2 = 0 and the

system spends no time in region V (in other words, system

does not stick, but is only momentarily motionless). As

discussed in Chapter III, the region V equations do not have

to be applied, since the entry and exit points for the

trajectory in this region are the same. If this also occurs

at the other switch point (i.e., T4 = 0), the limit cycle is

driven purely by sliding friction.

The other case is more complex to analyze, since the

trajectory must be calculated through all four switching

periods. If the driving torque is insufficient to overcome

the static friction at t = T1, the mechanical load stops for

a period T2 until the controller can integrate up and break

it loose. The system has come to rest in region III after

leaving region I, so system behavior is governed by the

model for region III during this period. Both static and

sliding friction effects will affect this limit cycle. For

this case, the system state evolves as

(4.7) x(t) = es (t-T1 (T) T1 < t < (T1 + T2)

from the initial point when the velocity went to zero at

t = T Note that no translation is necessary, since the

friction term is absorbed into the system model in the

static condition, represented by state matrix A

In either the pure sliding friction case (T2 = 0) or

the sticky case, the requirements on the system state at

t = T1 + T2 is that the driving torque is sufficient to

break loose:

(4.8) Yn(T1 + T2) = V'A x(T1 + T2) > LS (or = LS)

Since the sticky case involved the system evolving until

sufficient torque was developed to break away, the system

will move as soon as torque equals the static friction, so

(4.8) becomes an equality; note that the same comment

applies to condition (4.2b).

There is also the implicit requirement that the

velocity is zero when motion is impending (as in condition

(4.2)), but this requirement is not placed explicitly, since

it is taken care of by the model equations in region III

(which force acceleration to be zero) plus the zero velocity

at t = T1.

The system goes through the same behavior as above, in

the opposite sense, while traversing region II and returning

to III (or IV). The equations for the trajectory in these

regions, along with the conditions on the state at the

switching instants, are found by first performing the

translation to autonomous coordinates for region II, as in


(4.9) xla = x + (-Lc/an)


x = x + 1
-a -

i' = [ (-Lc/a ) 0 0 . 0 ]

as before. Note the friction torque is negative in this

region for LC > 0, since the velocity is positive. The

trajectory in region II is then

(4.10) x (t) = eA(t-T1-T2) x (T1 + T2)

= eA(t-T1T2 (x(T1 + T2) + i),

for T1 + T2 < t < T1 + T2 + T3

and the trajectory exits region II at the point where y

(velocity) again goes to zero

(4.11) x(T1 + T2 + T3) = x (T1 + T2 + T3) 1
= e 3 ( x(T1 + T2) + 1) 1

with the requirement that y (T + T2 + T) = 0:

(4.12) Yn(T1 + T2 + T3) = x(T + T2 + T3) = 0

At this point, the system trajectory either immediately

breaks free (trajectory in region IV) and leaves the region,

or it sticks (trajectory in region III) and follows the

trajectory defined by

(4.13) x(t) = e (t-T1-T2-T3) x(T1 + T2 + T3),

for T1 + T2 + T < t < T + T2 + T3 + T

from the point when the velocity went to zero. Note that no

translation is necessary, since the friction term is

absorbed into the system model in the static condition as in

the case of region III.

Collecting the conditions (4.5, 4.7, 4.11, 4.13) for

the trajectory through all four regions, and placing the

requirement that the end-point match the initial point, we

have the following nonlinear algebraic matrix equation for

the initial point:
(4.14) AT AT AT AT
(4.14) I es 4 e 3 eAs 2 e 1) x = 1

1 = es 4 (e 3 [e s 2 (I e A} 1 + 1] 1)
= -es4 e 3 eAsT2 e 1 1 + e s 4 e 3 e s 2 1
+ e s 4 e 3 1 es 4 1

In addition, there are four scalar nonlinear equations

from the four conditions on switching points (4.2a/b, 4.6,

4.8, and 4.12, but the last one is redundant). Altogether,

there are n+4 (nonlinear) equations in the n+4 unknowns x0'

T1, T2, T3, and T4. In the pure sliding friction case, of

course, T2 and T4 are both zero, but conditions (4.2b) and

(4.8) become inequalities and thus drop out, and (4.14)

simplifies to
(4.15) {I e 3 e 1) 0 = 11

1 = e 3 [ (I e 1 1 + 1 ] 1
= -eA3 eAT 1 1 + 2 eAT3 1 1

so there are n+2 equations in the n+2 unknowns xg, T1, and

T3 plus the two inequality conditions on the switch


The results presented so far in this chapter can be

summarized as

Theorem 4-1: A simple limit cycle exists for the system

model under study only if the following algebraic conditions

hold on the variables xO, T1, T2, T3, and T4:

Sticking Limit Cycle Case
(4.14) (I e s 4 e 3 e s 2 e 1) 0 = 1

where 1 = es 4 (e AT3 [es 2 (I e AT1) 1+ 1] 1)
= -e s 4 e 3 e s 2 e 1 1 + e s 4 e 3 e s 2 1
+ esT4 eA3 1 es4 1

(4.12) Yn(T1 + T2 + T3)

=[ 0 0 Pi n-3 1 ] x(T1 + T2 + T3)

= (T1 + T + T) = 0

(4.2b) Yn(0) = x' A x0 = -LS

(4.6) Yn(T1) = X(T1) = 0

(4.8) Yn(T1 + T2) = 2 'A x(T1 + T2) = LS

Pure Sliding Limit Cycle Case (T2 = T3 = 0)

(4.15) {I e 3 e 1) x0 = 1
-0 1
where 1 = e 3 [ {I e 1) 1 + 1 ] 1
= -eA3 eAT 1 1 + 2 e AT3 1 1

(4.2a) y (0) = [ 0 P . n-3 1 ] X

= X0 = 0

(4.6) Yn(T1) = x(T) = 0

In other words, the conditions above are necessary for the

existence of a simple limit cycle. Condition (4.12) was

included in the sticking case instead of the equivalent

condition (4.2a) to obtain a more symmetric form of the


Note that the conditions requiring velocity = 0 could

be applied at either t = T1 or t = T1 + T2 (equation (4.6)),

and at either t = 0 or t = T1 + T2 + T3 (equation (4.2a)).

This is a consequence of the fact that velocity remains zero

between these times.

Proof: The derivation presented above is a step-by-step

solution, through each region, for the limit'cycle

trajectory that was assumed to exist. Since the conditions

stated in the theorem were derived directly from the

trajectory solutions as seen above, and since a simple limit

cycle by definition must traverse those regions as assumed,

the conditions are necessary, and the proof is complete.

Note that a non-simple limit cycle, should any exist,

would not have to meet the necessary conditions. Thus, if a

solution to the necessary conditions exists for a system, a

limit cycle may exist; if no solution exists, no simple

limit cycle exists, but a non-simple one is not ruled out!

A complete theory requires an examination of the possibility

of more complex limit cycles.

It is interesting to note that the same development

could be performed using physical state variables. This

would result in replacing the control canonical A matrix in

the equations by the original system matrix, eliminating the

need for the vector P, and replacing 1 by the translation in

equation (3.12). The form used here clarifies the role of

system poles and zeroes, however, as seen by the following


The definition of the initial state X0, equation

(4.14), is completely determined by the poles (for fixed

switching periods), since the A matrix contains the

characteristic polynomial information only. Now, observing

that the transfer function from the friction input to the

velocity as an output is

Yn(s) / l(s) = [0 0 . 0 1] T12 (sI-A)-1 b
(this holds since y(s)=T12x(s), and x(s)=(sI-A)-lbl(s) ),

with b' = [0 0 . 1], it is clear that the transfer

function zeroes are defined by the bottom row of T12, that

is, the vector P. Thus, varying the zeroes only would

change the P vector only (in such a way that poles of the

closed loop system are unaffected), X0 would be unchanged,

and the four (or two) switching conditions alone would be

affected. The control form representation clarifies this

relationship, which suggests an interesting line for further


Illustrative Examples

Example IV-1: Two-Dimensional (2D) System with Sliding

Friction Limit Cycle

This is the same as Example I-1 of Chapter I; the

piecewise linear method is the same as standard phase plane

methods for this 2D case. The state variable model of the

system is

(4.16) dyl/dt = y2

dy2/dt = -(K/J) yl (B/J) Y2 Lf/J

where Lf is the friction torque, yl is position, and y2 is

velocity. The analysis performed below assumes K, J > 0,

although B can be negative (the other cases can be solved by

similar methods, but are left out for brevity). Chapter I

presents more information about the system, including

results of a previous solution by phase plane methods.

A second-order system cannot possibly have a limit

cycle with sticking, since once the system sticks, it is at

an equilibrium point (there is insufficient torque to break

free, and no controller integrator to ramp up). Therefore,

we need only examine the necessary conditions for the

sliding case. In addition, it will be assumed for this

example that the limit cycle is symmetric (this property is

proved for 2D systems in Chapter V). Applying the equations

developed in this chapter, the initial condition can be

found from equations (4.5) (which (4.15) simplifies to using

the symmetry assumption) and (4.6):
(4.17) x(T1) = eAT1 (x 1) + 1 = -

(4.18) y2(T1) = 1' x(T1) = x2(T1) = 0

where use is made of the fact that 0' = [ 0 1 ] (physical

state model already in canonical form), and x(T1) = -X0 (by


Making the assumption that K is large enough (or B is

small enough) that the system poles are complex, the matrix

exponential for this 2D system can be evaluated as

(4.19) eA = et cospt+(a/p)sinpt (1/p)sinpt

-(K/JP)sinpt cos3t-(a/p)sin3t

where a = B/(2J), 3 = (1/2) { (4K/J) (B/J)2 ).

Letting 0 = [a 0]' (since x2=0), and noting 1 = [LC/K 0]',

equations (4.17) and (4.18) become:
(4.20) x2(T1) = 0 = -[a (LC/K)] eT 1 (K/JP)sin3T1

which implies T1 = 7/3 (note a > (LC/K) since otherwise x0

is an equilibrium point).

Theoretically, T1 can also be any integer multiple of
7/3. However, noting that x2(t) = velocity = C e sin/t,

any multiple greater than one would entail sign changes in

the velocity for t
This is an example of a solution to the necessary conditions

that is not a valid limit cycle; see discussion on

sufficient conditions below.

The other component of x0 can be evaluated from (4.17)


(4.21) x1(T1) = -a
= [a (Lc/K)] eT l(cospT1+(a/P)sinPT1) + (LC/K)

= [a (LC/K)] e-a7/8 cosr + (LC/K)


(4.22) a + (LC/K) = [a (LC/K)] e-ar/

The initial condition can now be evaluated as

(4.23) a = (Lc/K) ((/+1)/(4-1)), P = e-ar/

which is the result presented in Chapter I. Since x0 is in

region IV, a must be greater than zero. This implies j>1,

requiring a<0, or B<0.

As the solution from Chapter I predicted, a sliding

type of limit cycle exists for every case with negative

damping (B<0) (if the position feedback gain K is

sufficiently large to give the linear system a complex pair

of poles). The initial conditions can be expressed in terms

of the system parameters by the equations from Chapter I:

(1.2) X10 = (LC/K) [(x+1)/(g-1)], x20 = 0


(1.3) p = exp[-Br/(2JP) ], = P [4K/J B2/J2]

Refer to Figure I-4 for a plot of amplitude (which equals

x10) as a function of gain K for the case J=LC=1, B=-l.

Example IV-2: Three-Dimensional (3D) System with Both

Sliding and Sticking Limit Cycles

The second example system to be considered (same as

Example I-2) is defined by the differential equations

(1.4) dy/dt = A y + b L


(1.5) A = 0 1 0

0 0 1

-K1 -K2 -B
with b' = [0 0 -1], Lf = friction torque, and where the

state variables are y2 = position, y3 = velocity, and

yl = compensator integrator. Figure I-6 shows a block

diagram of this system. The results presented in Chapter I

are now to be derived, based on the equations developed in

this chapter.

Case I: One Stable Pole and One Imaginary Pole Pair (Only

Sticking Limit Cycle Exists)

As in the example in Chapter I, let us first set

B = K = K = 1, and set the sliding friction torque to 1.0

and sticky friction (breakaway torque) to 1.2.

The equations (4.2a/b, 4.6, 4.8) defining the switching

periods and initial state (4.14) can be set up and solved

for this specific case, in order to demonstrate the

analytical calculation of the limit cycle trajectory. The

eigenvalues of this system are -1 and jl. The matrix

exponential of A is required for the equations, and can be

found (by, for example, diagonalizing A) as
AT -T -T
(4.24) 2eAT e +sinT+cosT 2sinT e +sinT-cosT
-T -T
-e -sinT+cosT 2cosT -e +sinT+cosT
-T s-T
e -sinT-cosT -2sinT e -sinT+cosT

The equations to be applied are
(4.25) x = e AT1 ( 1) + 1
(4.26) x2 = es 2 x -x0

where a symmetric limit cycle was assumed (breakaway point

at t = T1 + T2 is at -x0), and therefore the simpler

equations were used in place of (4.14). Note that A is the
A matrix with the bottom row set to zero, hence the matrix

exponential in the sticking region is

(4.27) eAsT2 = 1 T2 T22/2

0 1 T2

0 0 1

In addition, the breakaway condition that the acceleration

torque equal the static friction is required:

(4.28) 1'Ax0 = dx3(0)/dt = -LS

Noting that the translation vector 1' = [ 1 0 0 ], letting

0 = [ a b 0 ] (since x30 must be zero), and letting

xI = [ c d 0 ], we obtain

(4.29) c 1 = (1/2)(a-l)(e-T + sinT + cosT) + bsinT

d = (1/2)(a-1)(-e-T sinT + cosT) + bcosT

0 = (1/2)(a-1)(e- sinT cosT) -bsinT

from application of equation (4.25) and the definition of
the matrix exponential e,

(4.30) c + d T2 = -a

d = -b

from application of equation (4.26) and the definition of
e s 2, and

(4.31) a + b = L

from application of equation (4.28). By eliminating c and d

from these equations and performing additional algebra, they

simplify to

(4.32) b = LS a
T2 = [(a-l)e-T1 + (a+l)] / (LS a)
a [e 1 + sinT cosT ]
= eT1 + (21sf l)sinT1 cosT1
a [-e 1 sinT1 cosT1 -2]

= -e T1 (2 LS l)cosT1 sinT1 2LS

The first two equations define b and T2 in terms of a and

TI, while b and T2 were eliminated from the last two, which

can be used to simultaneously solve for a and T1.

Eliminating a and after some algebra, the last two equations

form a nonlinear algebraic equation, whose zeroes are

potential limit cycle solutions:

(4.33) f(T1) = (e-T 1 l)(l+cosT1) (e-T1 + l)sinT1 = 0

The plot of this function in Figure (IV-1) and analysis

shows that there are zeroes at odd multiples of r, and also

near 37/2 (+ 2n7, n = 0, 1, ...). The solution at r results

in a valid limit cycle solution where x = [ 1 0.2 0 ],
T1 = r, and T2 = 10, yielding an overall limit cycle period

of approx. 26.3 seconds, matching the simulation results.

The second potential solution is at T1 = 4.73 (approx. 37/2,

obtained by numerical solution of f(T1)=0), and is invalid,

since it results in a negative value for T2 (of -12.214).

Additional solutions are also invalid, for the same

reason as in Example IV-1. Examination of the matrix


__ _-------- ------- ,^.





O rr m

> 0 o 0 0





, -

Cn C

11C ti

exponential and the resulting formula for the velocity state

shows that the velocity changes sign approximately every 7

seconds. Therefore, the trajectory using these longer

values for T1 leaves region II before the switching time, so

it is invalid (see next section).

Therefore, the solution of the piecewise-linear

equations results in a valid limit cycle solution that

matches simulation results. Chapter I contains plots of the

limit cycle trajectory.

Case II: One Stable Pole and One Unstable Pole Pair (Both

Sticking and Sliding Limit Cycles Exist)

If the same example system is used with B = 0.9, and

K1 = K2 = 1 still (case II in Chapter I), the closed-loop

eigenvalues (of the linear portion of the system) move into

the right-half plane, to .026 j 1.024, while the other

pole is at -.9524. As stated in Chapter I, there is still a

sticking limit cycle solution close to that of Case I

(Figures I-7 through 1-9, confirmed by simulation). A

similar analysis to that presented above can'be used for

this case, so the details are omitted.

Note, however, that numerical methods may be required

to evaluate the matrix exponential (due to the tediousness

of the calculations; case I was a particularly simple form);

this forces the use of iterative methods to calculate the

limit cycle parameters (see example of this method below).

The second limit cycle solution, as discussed in

Chapter I, is of the sliding type. Simulation having shown

one with a half-period T1 approximately equal to 3.1

seconds, an iterative numerical solution of the equations

for the symmetric, sliding case was performed, as follows.

The equation for this case is obtained from (4.5),

where again the assumption is made that the limit cycle is

(4.34) xl = eAT ( 0 1) + 1

= -xO

which can be solved as

(4.35) xO = (I + eAT1)-1 (I eAT) 1

A computer program calculated values for the initial

condition xg, given trial values of T1, and the process was

repeated until the velocity initial condition XO(3) = 0.

This was not difficult, since an approximate starting

condition was available from the simulation results. Given

a T1 that resulted in a zero initial velocity, the entire

initial condition could be defined from (4.35).

This yielded a solution at T1 = 3.1455 seconds, and

20 = yO = [-0.283 12.4 0.]. The solution is valid since

the breakaway condition is also satisfied (torque at zero is

greater than LS). Figures 1-10 and I-11 show some views of

the limit cycle for this case. An appendix contains

computer code that was used to perform this iterative

calculation. The code applies to the sliding case only, and

would have to be modified for the general case.

Exact (Necessary and Sufficient) Conditions

Theorem 4-1 presented some conditions that were

necessary for the existence of a simple limit cycle. That

they are not sufficient is seen by the results in the

Examples, where solutions were found that met all of the

necessary conditions, yet were not valid limit cycles (the

switching period T2 was less than zero, or velocity changed

sign during a switching period).

In order to expand the previous conditions into a set

that is also sufficient, the concept of a consistent

solution is useful:

Definition: A consistent solution is a solution of the

equations in Theorem 4-1 that meets the assumptions about

the region containing the trajectory during each switching

period. These assumptions are

(Sticking Case)

(1) Trajectory in region I, 0 <

(2) Trajectory in region III, T

(3) Trajectory in region II,

T1 + T2 < t < T1 + T2

(4) Trajectory in region III,

T1 + T2 + T3 5 t T1

(Pure Sliding Case)

(1) Trajectory in region I, 0 <

(2) Trajectory in region V, t =

(3) Trajectory in region II, T,

t < T1

1 t 2 T1 + T2

+ T3

+ T2 + T3 + T4

< t < T1 + T3

(4) Trajectory in region IV, t = T1 + T3

A check of the consistency of a solution is

straightforward, based on the definitions of the various

regions presented in previous chapters: the sign of the

velocity is examined for assumptions (1) and (3), and the

magnitude of the torque for (2) and (4). Given this

definition, a set of necessary and sufficient existence

conditions can now be presented:

Theorem 4-2: A simple limit cycle exists in the system

under study if and only if a consistent solution of the

Theorem 4-1 equations exists. In addition, this limit cycle

has the exact properties (switching periods, initial

condition, and trajectory) defined by the corresponding

solution to the equations.

Proof: (Necessity) Suppose a simple limit cycle exists

with the given properties. Theorem 4-1 (more exactly, the

derivation in the first part of this chapter) demonstrated

by analysis of this limit cycle that the listed equations in

fact apply to the trajectory.

In addition, this derivation, along with the assumption

that the limit cycle is simple, shows that the solution must

meet the consistency conditions, i.e., it must be in the

appropriate region of the state space during each switching

period. Therefore this condition is also necessary.

(Sufficiency) Now suppose a consistent solution of the

equations exists. A limit cycle can be generated, with the

same properties as defined in the solution, by starting at

x0 and following around to the initial point again. At each

step of the process, the consistency conditions require the

trajectory to be within a given region, so that the system

equations for that region apply. The equations listed in

Theorem 4-1, derived from the application of those system

equations, then show that the trajectory reaches the next

switching point at the appropriate time and state.

Therefore, a limit cycle must exist with those exact

properties, and the proof is complete.

Although this theorem defines a set of exact

conditions, i.e., they are equivalent (necessary and

sufficient) to the existence of the limit cycle, it is

appropriate to ask if these conditions are minimal. In

other words, are there redundancies in the set of

conditions? Is there a simpler set of conditions that are

still exact? For example, it may be sufficient to merely

require switch periods greater than zero, thus eliminating

cases such as found in Example IV-2. These questions are

open, and may be subjects of further research.


This chapter uses the results of chapter IV and other

known facts about the friction nonlinearity to examine the

behavior and characteristics of friction limit cycles. The

model of the system developed previously is used to analyze

the stability of limit cycles predicted by the equations in

chapter IV. In addition, the relationships between the

oddness of the system differential equations, symmetry, and

uniqueness of limit cycles is explored.

Stability of Predicted Limit Cycles

Once a limit cycle has been predicted by the solution

of the equations derived in Chapter IV, it is natural to ask

about its relationship to the global phase portrait. Since

the piecewise linear model used here is exact, the exact

limit cycle period, amplitude, and trajectory can be found

immediately once an initial point on the limit cycle is

known. It is more difficult to determine stability and

other information about the phase portrait, however, since

it can require extensive simulation.

However, local stability of the limit cycle can be

determined by the standard method of linearization of the

periodic flow map at the initial point previously found. In

other words, an approximate model of the trajectories near

the limit cycle initial point can be determined by


linearization. This analysis provides information about the

asymptotic stability of a closed orbit solution as discussed

in Hirsch and Smale (1974), chapter 13 (if more information

on orbital stability analysis is of interest, this reference

provides a good description, and was quite useful in the

development of the results in this chapter).

Note that the model provided by this analysis is not a

continuous-time representation of trajectories near the

cycle; instead, this discrete-time model shows the deviation

from the initial point and how it evolves after each cycle

period. Specifically, an initial sufficiently small

deviation from the initial point x0 is assumed (so the

trajectory starts at x +6x ) and the trajectory traced once

around the limit cycle. The new deviation from the initial

point when the trajectory completes the cycle is related to

the initial deviation by the equation

(5.1) 6xT = Y 6Sx + higher order terms

where 6xT is the new deviation from the initial point

(0+6xT is the point on the trajectory after one cycle), and
Y is the transition matrix, representing the first-order

evolution of the deviations from the initial point over

time. Eigenvalues of the transition matrix determine

stability, while eigenvectors determine the phase portrait

near the limit cycle.

This analysis is attempted below; however, because of

the discontinuities in the system nonlinearity the analysis

is difficult. Therefore, the analysis is accomplished, as

in the derivation of Chapter IV, by solving one region at a

time, and then pasting the trajectories together. The first

case analyzed is that where the predicted limit cycle has no

sticking (sliding friction only). The case with sticking

during the cycle is then analyzed by building on the first


Case 1: Limit Cycle with no Sticking

Assume an initial point x0 that meets the equations in

Chapter IV with no sticking. In this case, the first

switching period is equal to the constant T1 in chapter IV,

the second switching period is T3, and the constants T2 and

T4, representing the periods of stickiness, are zero.

A deviation 6x0 is assumed from the initial point of

the limit cycle. We wish to limit deviations, however, so

that the perturbed initial point is still contained in

region IV. This is necessary in order to apply the Chapter

IV equations, since the initial point to which the equations

apply is assumed to be in region IV. Limiting the deviation

in this fashion results in no loss of information about the

limit cycle stability, however; (n-l) of the eigenvalues of

the periodic flow map can still be determined, while the nt

eigenvalue (that applies to deviations out of region IV, and

along the limit cycle) would be zero, as will be seen in the

examples below. Note that the derivative of the flow after

one orbital period has an eigenvalue of unity (Hirsch and

Smale, 1974, p. 277), but the Poincare map has a zero

eigenvalue (Hirsch and Smale, chapter 13, section 3).

In order to limit the deviation 6x such that the

initial point remains in region IV, 5x can be constructed

-1 '
(5.2) 6x = T12 Ir 6y0r

where T12 is the composite linear transformation from

Chapter III used to convert between canonical and physical

state variables (y = T12 x), I is an n x (n-l) matrix (Ir

is (n-1) x n) used to force the nth element of the deviation

vector Sy0 to be zero:

(5.3) I = 1 0 0 ... 0

0 1 0 ... 0

0 ... 0 1 0

0 ... 0 0 1

0 ... 0 0 0

and, finally, Sy0r is the (n-l)-vector of deviations in

physical state variables. This construction forces the 6x0

vector to remain in region IV, since the change in initial

velocity (6Sy(n), the nth element of the physical state

vector) is zero. Although the derivation below will be done

in canonical variables, this construction (5.2) will be used

at the end to complete the calculation for physical

variables, and force the deviation to have zero initial


Starting at the initial point x +6x the trajectory is

determined by the piecewise linear methods as in chapter IV.

The trajectory in region I is found using autonomous state


(5.4) x (t) = x(t) 1

so the trajectory is
At At
(5.5) x (t) = e x (0) = e [x0 + x-
-a -a
x(t) = e [X0 + SX0 1] + 1

For a sufficiently small deviation 6x0, the trajectory

reaches region V (velocity goes to zero preparatory to

reversing direction), but not necessarily in exactly T1

seconds. Denoting the change in this period is ST, the

velocity goes to zero when the trajectory reaches

(5.6) x(T1 + 5T) = eA(T1+6T) [x + 5x 1 1

The quantity of interest is the deviation of this point from

the point where the original limit cycle reaches region V
(5.7) x e 1 [X 1] + 1

(5.8) Sx = x(T1+6T) x-
=(e I} e 1 [ ] + e(T1T) 6x

Note that the first term shows the dependence on ST, since

it is zero if ST is zero, while the second term shows the

dependence on S6x0
Now that the exact value for Sx is determined, this

nonlinear formula must be linearized to determine local

behavior. One factor that causes difficulties is that 6T

depends on 6x0, which complicates the calculation for

(5.9) Y = (d(Sx )/d(5x )}|Ix=

where Y is the desired transition matrix defined above, and

the derivative is evaluated for 6x = 0. The dependence of

Sx on both ST and 6x0 explicitly is handled using the chain


Note that since x is not a function of 6x or ST,

(5.10) [d/d(6x0)] (6x ) = [d/d(6x )] {x(T1+6T) x)

= [d/d(Sx0)] (x(T1+6T))

= [d/d(6x0)] (eA(T1+6T) [X0 + S 1])

where the term in x(T1+6T) equal to 1 was also dropped since

it does not depend on x 0. Using the chain rule, this

derivative is

(5.11) d(6x1)/d(6x0) = D(6X1)/D(6x0)

+ [D(6X1)/D(6T)] [d(6T)/d(6x0)]

where the symbol "D" indicates a partial derivative. The

first term is

(5.12) D(6x )/D(6x0) = eA(T1+6T)
= eAT1

when evaluated at 6x0=0 (hence 6T=0). The partial

derivative in the second term can be evaluated as follows.

First define
m = e 1 [X0 + &x 1]

Then, by using (5.10), the desired partial derivative is

(5.13) D(6x1)/D(6T) = [D/D(6T)] {eAST M)
A e m
= Ae 1 [ 1]

when the equation for m is substituted back in, and the

derivative is evaluated at 6x0=0 (6T=0). The

differentiation step can be verified by expanding the

exponential into a power series, and evaluating term by


The evaluation of the [d(6T)/d(6x0)] factor in the

second term requires the definition of a functional

relationship between 6T and 6x0. We know from the

derivation of x(T1 + 6T) above that the velocity goes to

zero at this point (trajectory reaches region V). As in the

analysis in Chapter IV, this can be stated as

(5.14) P' x(T1 + ST) = 0

(eA(T +6T)
= {e 1 [x + 6x0 1])

since P'1 = 0. This formula implicitly defines the

relationship between 6T and 6x in that we have a function

of the form f(6T,6x0) = 0. The implicit function theorem

(Rudin (1976), pp. 223-8) can be used to evaluate the

desired derivative as
(5.15) d(6T)/d(6x0) = [Df/D(6T)]-1 [Df/D(6x0)]

Evaluating these factors,

(5.16) Df/D(ST) = [D/D(6T)] {('(eA(T1+T) [0+6x0-11))
= [D/D(6T)] {(' eAST m)
= A e A

with m defined as above (to see this, again expand

exponential into series and evaluate term by term). The

second factor is

(5.17) Df/D(6x0) = [D/D(6x0)] {('(eA(T1+6T) [x0+60-]))

= [D/D(6x0)] {(' eA(T1+6T) 6x0)

= eA(T1+6T)

The derivative from the implicit function theorem, (5.15),

can now be evaluated as
AT -1 AT
(5.18) d(6T)/d(6x0) = -({'A eAT1 [O-I]} ) ('eAT 1

where the factors are again evaluated at 6x =0 (6T=0).

Finally, combining (5.12,13,18) into (5.11), we get the

total derivative representing the linearized periodic map

(5.19) Y = (d(6xl)/d(6x0)) 16x=o
= eA1 (A eAT1 ['X-1]) ('A eAT1 [-O i-{e A1)

Note that the factor to be inverted is a scalar. This

expression can be simplified by noting that
(5.20) i = eAT1 [x 1] + 1

so substituting into (5.19) gives
AT -1 AT
(5.21) Y = eAT1 (A [X1-1]) (P'A [X1- i)-1 ('eAT1)

= eAT1 (A [xl-1] 'e AT1) / ({'A [x-1])

The linearized flow map for the entire cycle must now

be constructed, using this equation (5.19) for the

linearized flow map for a half-cycle. Note, first of all,

that the linear transformation representing the linearized

flow map for the cycle is the product of the'linear

transformations of the two halves. This is a consequence of

the chain rule, since

(5.22) 6xT = f (6x ) = (f2(6x0))

d(6xT)/d(6xo) = (d(6)/d(6x1)) (d(6x1)/d(6x ))

The second of these two linear transformations has already

been derived; the first must be found and the linear map for

the entire cycle found by multiplication.

Assume first of all that the limit cycle is symmetric,

so that the point at which the limit cycle velocity goes to

zero is exactly one-half period after the initial point, at

-X0. In this case, the half-period is equal to the constant

T1 and T3 equals T l

Therefore, the linearization in which we are interested

is completely determined by the first half of the limit

cycle trajectory, because of the odd symmetry. To show

this, assume that the same deviation from the point xl is

taken, except with a change of sign, and trace the

trajectory from this point x -6x Noting that this initial

point equals (-x0-6x0), and using the odd symmetry of the

phase portrait, it is clear that the trajectory returns to

region IV at the point (-x -6x ) (since the trajectory from

x0+6x0 goes to x1+6x ). Thus the resulting deviation on the

second half of the cycle is

(5.23) 6x = (-x2-6x ) X

= -6X1

since xl = -x0 Linearizing this second leg, using the

initial deviation of -6x ,

(5.24) 6x = Z (-6x ) + higher order terms
= -6x

= -Y 6x + h.o.t.

by the analysis of the first half-cycle. Since this holds

true for any deviation 6x0, we must have Y = Z. Therefore,

the transition matrix for the whole cycle is the square of

the transition matrix for the first half-cycle, Y2

Now remove the assumption of symmetry. Assume that the

limit cycle is not symmetric (but still involves only

sliding, i.e. case 1, no sticking); in this case the

linearized flow map must be examined on the second half of

the cycle also. The same derivation as above (equations

(5.4) through (5.8)) can be performed, but starting at the

initial point x1+6x1, in order to determine SxT in terms of

6x1. The trajectory in region II is found using autonomous

state variables:

(5.25) x (t) = x(t) + 1
so the trajectory is

(5.26) x(t) = eA(t-T) [ + x1 + 1] -1

For a sufficiently small deviation 6x1, the velocity goes to

zero when the trajectory reaches

(5.27) x(T + 6T3) = eA(T3+6T3) [x + 6x + 1 -

The deviation of this point from the point where the

original limit cycle reaches region IV is

(5.28) 6xT = x(T+ST3) xT

= (eAT3 I) eA3 [x + 1] + eA(T3+T3x

Note that by comparing this formula to (5.8), and stating

(5.8) as a function with parameters

(5.29) 6x= (eAST I) eAT1 [x 1] + eA(T1 T)

= f(6x0, 6T; T1, x0-1)

then we see

(5.30) 6x = f(6xl, ST3; T, x +1)

The linearization for this function has already been

derived, so the linear map for the second half of the cycle

is found by substituting the correct parameters from (5.30)

into (5.19), so that

(5.31) 6ST = Y(T3, x1+1) 6x, + h. o. t.

(5.32) Y(T3, l+1) = d(6XT)/d(6X1) 6x=0
= eAT3 (A eAT3 [+]) (~'A eAT3 [+)-1 ('eAT3)

By the argument given previously in (5.22), the total

linearized flow map is the product

(5.33) Ycompos = Y(T3' +1) Y(T1' 0-)

6T = Y cx + h. o. t.
-T compos -0
Note that this reduces to the symmetric case when T3 = T1

and x = 0:

(5.34) Ycompos Y(T -x+1) Y(T, xO-1)
= Y2(T1, -1)

since Y is an even function of its second variable.

The final step in the derivation for this case is to

convert back to physical variables, while at the same time

forcing the initial deviation to the desired region. By

following the reasoning given at the start of the analysis

for this case (equation (5.2)), the desired matrix for the

linearized flow map of the entire cycle is

(5.35) Y -compos = Y(T' +1) Y -)
r-compos r 3' l1 O ( 0
Yr(T, x -1) = Ir T 2 Y(T x0-1) T Ir

(5.36) 6YTr = Yr-compos 60r + higher order terms

where 6yTr is the (n-l)-vector representing the deviation in

physical variables after going around one cycle (with

velocity = 6y left out since it is zero), 6y0r is the

initial deviation in the same variables (also has (n-l)

elements, no deviation in initial velocity allowed), and Yr'

Y are the matrices representing the linearized flow
maps derived above, but in physical variables and reduced to

order (n-1). Note that the latter matrix is (n-l)x(n-1), so

its (n-1) eigenvalues are those desired for the local

stability analysis (as discussed above, the nth eigenvalue

is always 0).

Case 2: Limit Cycle with Sticking

To analyze this case, it is necessary to derive the

linearized flow map for the portions of the limit cycle in

the sticky region, i.e., region III. These linear maps can

then be combined with the maps previously derived for the

sliding regions to form the composite map for the cycle.

Assume an initial point x3 that meets the equations in

Chapter IV, so that this is the point at which the system

comes to rest in region III (this is the fourth switch point

in the cycle). In this case, the period the system sticks

in region III is equal to the constant T4 in Chapter IV.

This initial point must be used, instead of the point x at

the point of impending motion, in order to derive the

behavior within region III.

A deviation 6x3 is assumed from the initial point of

the limit cycle. This deviation is of course constrained as

before, so that the resulting initial point remains inside

region III (velocity = 0). Starting at the initial point

x3+6x3, the trajectory is determined by the piecewise linear

methods as in chapter IV.

The trajectory in region III is

(5.37) x(t) = eAt x(0) = eAt [X3 + 6x 3

For a sufficiently small deviation 6x3, the trajectory

reaches the condition of impending motion (acceleration

torque sufficient to break loose), but not necessarily in

exactly T4 seconds. Denoting the change in this period is

ST, impending motion is reached when

(5.38) x(T4 + ST) = eAs(T4+6T) [x3 + 6x3]

The quantity of interest is the deviation of this point from

the point where the original limit cycle breaks free
(5.39) x es 4 x

(5.40) 6x = x(T4+6T) 4
A 6T AT A (T +6T)
= {es I) es 4 x + esT4T)
The differentiation of this function to determine the

first-order dependence of 6x on 6x3 is performed as before;

many of the steps are similar and thus omitted. The desired

derivative is

(5.41) Ys = d(6x4)/d(6x3) 1=0

where Y is the desired linearized map in the stiction

region. Note that since x4 is not a function of 6x3 or 6T,

(5.42) [d/d(6x3)] (6x4} = [d/d(6x3)] {x(T4+6T) x4)

= [d/d(6x3)] (x(T4+6T)}

= [d/d(6x3)] {eAs(T4+6T) [x3 + 6X3])

The dependence of 6x on both ST and 6x3 explicitly is

handled using the chain rule, as before:

(5.43) d(6x4)/d(6x3) = D(6x4)/D(6x3)

+ [D(6x4)/D(6T)] [d(6T)/d(6x3)]

where "D" is again the symbol for partial derivative. The

first term is

(5.44) D(6x4)/D(6x3) = e s

when evaluated at 6x = 6T = 0. The partial derivative in

the second term is
A 6T
(5.45) D(6x4)/D(6T) = [D/D(6T)] {e s m)
= A e s 4 [x3]

where m is a temporary vector independent of 6T, as in the

previous derivation.

The evaluation of the [d(6T)/d(6x3)] factor in the

second term is handled as before by relating the two

variables with an implicit function, and using the implicit

function theorem. The implicit relation in this case is

determined by the breakaway condition (i.e., impending

motion); as in the analysis in Chapter IV, this can be

stated as

(5.46) P' A x(T4+6T) = -LS

or, equivalently,
0' A eAs(T4+6T) [x3 + 6x3] + LS = 0

Thus we have a function of the form f(6T, 6x ) = 0. The

derivative is then

(5.47) d(6T)/d(63) = [Df/D(6T)]-1 [Df/D(6X3)]

Evaluating these factors,

(5.48) Df/D(6T) = [D/D(6T)] ({' A e s6T m)
A 6T
= e' AA es m

The second factor is

(5.49) Df/D(6x3) = [D/D(6x3)] ({' A eAs(T4+6T) 6 x

= A eAs(T4+6T)

so that
AT -1 AT
(5.50) d(6T)/d(6x ) = -{('AA eAs 4 x3) {('A eAs 4)

where the factors are again evaluated at 6x =6T=0. Finally,

combining (5.44,45,50) into (5.43), we get the total

derivative representing the linearized periodic map

(5.51) Ys(T4, x3) = d(6x4)/d(6x3) 6x=0
= e s 4 (As e s 4 x3) ('AAs e s 4 x ) ('A e s 4)

Note that the factor to be inverted is a scalar. This

expression can be simplified by noting that
(5.52) x4 = e s x3

so substituting into (5.51) gives

(5.53) Ys(T4x3) = eAsT4 (As x4)('AAs 4)('A e s4)
= e s 4 (As x4 'A e s 4) / {('A A x4)

In order to obtain the linear map for the entire cycle

for the case with sticking during the cycle, we can use the

chain rule as before:

(5.54) Ycompos = d( )d(/d(6 2)) (d(6x2)/d(6x ))

x (d(6x )/d(6x0 ) (d(6x4)/d(x3 )

where the last factor was just derived, the first and third

are obtained from the sliding case previously analyzed, and

the second factor can be found from the fourth by symmetry,

as shown below. Note that 6x0 = 6x, since both are the

deviations at the point of impending motion in region III;

this equivalence justifies the juxtaposition of the last two

factors when the chain rule is applied.

Note also that in (5.54) the order of the factors is

different from that used in the sliding friction case; the

initial deviation for the case with sticking is taken as 6x3

at the point x where the system comes to rest, rather than

at the point x where motion is impending. This order is

required since any (small) deviations Sx3 are allowed that

keep the initial velocity zero (hence 6x3 can be constructed

as in (5.2)). On the other hand, certain initial

perturbations would not be allowed at the point x ; any

perturbations at this point of impending motion that altered

the forcing torque would either cause motion (so initial

point no longer in region III), or would mean motion is no

longer impending. This problem does not arise in the

sliding friction case, but in the sticking case there are

two constraints on x0: zero velocity and forcing torque at

breakaway. It is therefore much more convenient to take the

initial point to be x3 for the stability analysis; all

constraints are then taken care of in the derivation


The second factor in (5.54) can be found by symmetry,

as follows: taking a deviation from the point x~ at which

the system comes to rest in region III of Sxl, we can derive

by the same steps as above that

(5.55) 6x2 2=(T2+6T) x

=(es I) eAsT2 x + eAs(T2+T)65
1 1

The similarity of this equation to (5.40) can be used to

show that the desired linearization is

(5.56) Ys(T2, 1) = d(6x2)/d(6 ) 6X=
= es 2 (As es 2 xl} ('AAs e s 2 xl)1 'A eAsT2}

Of course, as in the sliding case, this matrix must be

modified to transform to physical variables, and reduced to

order (n-l).

The linear transformation for the entire, four-part

cycle is then

(5.57) Yr-compos = Y(T3' 2+1) Y (T x )
r-compos r 2 rs 2'1
x Yr (T1,X-i) Yrs(T4'x3)
Y (T1' 0x-1) = I T2 Y(T -1) T2 I
12 r Y1' -0 12 I
Y(T2' l) =Ir T 2 Ys(T2' l) T12 I

(5.58) STr = Yr-compos XOr + higher order terms

where, as before, the Syr vectors are the (n-l) elements of

the initial and final deviation, respectively, and Y
is the desired linear transformation that determines local

stability of the limit cycle.

The results of the stability analysis of this chapter

can be summarized as follows:

Theorem 5-1: Given the existence of a friction limit cycle,

as predicted by the solution of the nonlinear algebraic

equations derived in Chapter IV (Theorem 4-1), the local

asymptotic stability of the periodic orbit is determined by

the eigenvalues of the matrix:

(Case 1: Sliding only)

(5.35) Y compos = Y (T +) Yr(T -)
r-compos r 3' -1 r(1' -0

(Case 2: Sliding plus sticking)

(5.57) Yr-compos = Y(T3' x2+1) Y rs(T2, 1)
r-compos r 3' 2 rs2' x1
x Y (T1, xO-1) Yrs(T, x3)

Yr(T XO-1) = I1 T12 Y(TI, X--1) T12 I1
Yrs(T x) = I T2 Y (T, x ) T2- Ir
the x 's are the switching points, the T 's are the

corresponding switching periods, Ir is the (n-l) x n matrix

defined in (5.3), T12 is the transformation from physical to

canonical variables defined in Chapter III, and the matrix

functions Y and Y are
(5.19) Y(T1, x0-1) = d(6Xl)/d(6x0) 6x=0

SeAT1 (A eAT1 [X0-i]) ('A eAT1 [x -1])-{('eAT1)

(5.56) Ys(T2, X1) = d(6x2)/d(6x1) 6x=0
= es 2 (As e s 2 X1) ({'AA eAs 2 X}) ({'A es 2

Before this theorem is proved, some terms must be

defined more rigorously.

Definition: The flow of the differential equation system is

the map (t,x0):'(t,x0), where 9(t,x0) = x(t)'when the

initial condition x(0) = x .
Definition: The Poincare map (P-map) is the map 6x :g(6x ),

where g(6x0) = x(t) x0, 6x0 = x(0) x, and tl is the

first time at which the trajectory again crosses the local

section around the periodic orbit at x .
Proof: The proof follows from the theory in Chapter 13 of

Hirsch and Smale (1974) and the calculations above, except

for one difficulty to be overcome. The reference theory is

stated for the case where f(.) is continuously

differentiable (so that the flow < is continuous). However,

as will be shown below, all that is required to obtain the

desired result is that the flow be a continuous function of

the deviation about the initial point, x + Sx so that

certain limits can be taken, and, in addition, that the

P-map is differentiable, so that the derivative may be used

to check asymptotic stability.

As a preliminary step, then, these properties must be

shown. The flow is the integral of a real-valued function

of time (see equation 5.69 below), f(x(.)), so it is a

continuous function of time (Rudin, (1976), Theorem 6.20).

It is also clear that the flow is continuous in x within

each piecewise linear region, since the function f(.) is at

least C1 there (and Hirsch and Smale, Chapter 15, shows this

is sufficient for continuity of the flow). Now, taking the

initial condition as a sufficiently small deviation from the

initial point x0' the flow at any finite time is the

composition of a finite number of continuous functions,

formed from the pieces of the trajectory in each region it

traverses. This composition can be formed since the

trajectory is continuous across the discontinuities in f(.);

that is, x(t+) = x(t-), where t is a switching time.

Therefore the flow is a continuous function of this initial

condition (Rudin, (1976), Theorem 4.7, the composition of

continuous functions is continuous), for sufficiently small


The same approach is used to demonstrate the

differentiability of the P-map. This map is simply the

composition of the flows in each of the regions, so is

continuous by the previous argument. But, in fact, the

flows in each region are continuously differentiable, so

their composition is also differentiable.

Using these properties, the rest of the proof proceeds

as in Hirsch and Smale (1974). In chapter 13, section 3, it

is shown that if the P-map has x0 as a sink, then the orbit

is asymptotically stable (using the property of continuity

of the flow to show that the flow converges uniformly to the

periodic orbit if the deviations after each cycle also

converge). But the P-map is a discrete-time dynamical

system, so it is a sink if the derivative of this map has

all eigenvalues less than one in magnitude (which requires

the property of differentiability of the P-map). The

calculations of this derivative have been presented above,

and show that the matrix Y is in fact the derivative
of the P-map at x Therefore the (n-l) eigenvalues of this

matrix determine if the orbit is asymptotically stable, and

the proof is complete.

To summarize the above stability results, it has been

shown that the perturbations from the periodic orbit evolve

according to the equation

(5.58) yTr = Yr-compos 65Or + higher order terms

where 6YTr is the (n-l)-vector representing the deviation in

physical variables after going around one cycle (with

velocity = Syn left out since it is zero), and Sy0r is the

initial deviation in the same variables (also has (n-l)

elements, no deviation in initial velocity allowed). The

eigenvalues of the matrix Yr-compos' therefore, indicate the

stability of the periodic orbit (i.e., if they are all less

than one in magnitude, the orbit is asymptotically stable,

and in fact a periodic attractor).

In addition, it has been shown that the results on

orbital stability in Hirsch and Smale (1974), chapter 13,

can be extended to the case where f(.) is only piecewise

linear, rather than C1.

Stability Calculations for Example Problems

The limit cycle orbital stability equations derived

above can be applied to the example system discussed in

Chapters I and IV (example 1-2, case 2, a 3D system with

both a sliding and stiction limit cycle). The example

system is defined, as in the example, by

(5.59) dy/dt = A y + b Lf


A = 0 1 0

0 0 1

-K1 -K2 -B
with b' = [0 0 -1], Lf = friction torque, and where the

state variables are y2 = position, y3 = velocity, and

y = compensator integrator. The transformation matrices T1
and T2 are identity matrices, and the row vector 9' equals

[ 0 1]. The particular case has B = 0.9, K = K2 = 1, with

closed-loop poles at .026 j 1.024 and -.9524. The

stiction limit cycle solution is at approximately

x = [0.98 0.22 0], with T1 at approximately 7. There is

also a solution of the equations for the symmetric, sliding

case, at T1 = 3.1455 seconds, and xO = [-0.283 12.4 0].

Example V-l, Part 1: Stability of Sliding Limit Cycle

For this case, equations (5.35, 5.19) from Theorem 5-1

apply. Using the parameters defined above, we find that
(5.60) Ae 1 [x 1] = [-12.41 .001 13.13]'

so that

P' AeAT1 [x0

The matrix exponential

algorithm was used)

- 1] = 13.13

is (approximately,

(5.61) eAT = -.5314



Therefore, the matrix Y required

flow map is




to define

a numerical




the linearized

(5.62) Y = .0197 -.0556 .0223

-.5227 -1.054 -.583

.0002 .0000 -.0002

Although there are unavoidable numerical errors in this

calculation, the true Y matrix apparently does have a zero

eigenvalue, as expected. The upper left block sub-matrix is

equal to Yr' and the composite stability matrix is

(5.63) Y = .02944 .05749

.5407 1.14

This matrix has eigenvalues of .002125 (with corresponding

eigenvector at [-2.104 1.]' ), and 1.168 (with eigenvector

[.05 1.]' ). Thus it has one stable and one unstable mode

(note this is a discrete-time system), and is a saddle

point. A simulation of this system was used to verify this,

where the xl and x2 states were sampled whenever the

velocity x3 went through zero (near x0). The plot in Figure

V-1 shows the point x0, the eigenvectors of the stable and

unstable modes, and typical trajectories near x The

stable mode is very fast so that simulated trajectories

quickly jump onto the unstable eigenvector, then move along

it. The simulation experimentally verifies the stability


Example V-l, Part 2: Stability of Sticking Limit Cycle

Similar calculations can be performed for the other

limit cycle that exists for this system, except that

equation (5.57) from Theorem 5-1 governs this case. Using

the initial condition defined above, with T1 = v and T 9

(as approximate values), and xl = x(T1) = [1. -0.22 0.]',

and recalling that the symmetry defines the other switching

points and periods, we find after substitution in (5.56)

(5.64) Y = 0 -1 -9.9

0 1 9

0 0 1

with Y being the upper left block sub-matrix. The

stability matrix for the sliding portions of the limit cycle

is also required, and the upper left block of Y is

University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - - mvs