“mom.

In";

.3031 ,
cu v

y uVJw 1.x

..w .x .‘r.

J yx.'z
31:34,

.
t
a

1: 1.1...
571.125.:

 

In...
hat...

:5

Wu!

 

 

MRI-n...

:Ibit'éuu

nun

z¢iV

mum:

I!’
41‘

a! a. .fhuii.

uhhuf

 

THESIS WNW““1111111“WW‘MHIHIUUI\IW

:3 3 1293 014172

7“»

This is to certify that the
dissertation entitled
A METHOD OF PARAMETRIC IDENTIFICATION
FOR CHAOTIC SYSTEMS

presented by

Ching-Ming Yuan

has been accepted towards fulfillment
of the requirements for

Ph.D. , Mechanical Engineering
degree in

 

 

M.

r

 

dajor professoj

Date ‘61 IN}, 9(

MS U is an Affirmariw Action/Equal Opportunity Institution 0- 12771

 

way wa— " m:

 

 

LIBRARY
Mlchlgan State
Unlverslty

 

 

 

PLACE N RETURN BOX to remove this Moot from your record.
TO AVOID FINES Mum on or before date duo.

DATE DUE DATE DUE DATE DUE

       
   

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

I’T’I s

MSU In An Affirmative Action/Equal Opportunity Institution
W

A METHOD OF
PARAMETRIC IDENTIFICATION FOR
CHAOTIC SYSTEMS
By

Ching-Ming Yuan

A DISSERTATION

Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of

DOCI' OR OF PHILOSOPHY

Department of Mechanical Engineering

1995

ABSTRACT

A METHOD OF
PARAMETRIC IDENTIFICATION FOR
CHAOTIC SYSTEMS

By

Citing-Ming Yuan

Amethodforidenfifyingpmmeminamathemadealmodelofachaodcsystemis
presented. It is an extension of an existing method for nonlinear systems with Stable
periodic response. The method exploits the chaotic attractor, and extracts the unstable
puiodicmbitsnomtheamactormteprecentthesystembehavim.fiachtaminthe
mathematical model is expressed in a finite Fourier series using the extracted periodic-
orbits,andthe harmonic-balancemethodisappliedtoformasetoflinearalgebraic
equations in system parameters for least-squares estimation.

This method has been successfully applied numerically to a forced Duffing oscillator, a
smooth Coulomb friction system. a parametrically forced system. and a Lorenz oscillator,
and experimentally to a forced oscillator with a two-well stifi'ness potential.

The identified models have been verified by comparing the Lyapunov exponents. the
suucuneoftheummblepaiodicorbimandmebifincadondiagramsofdwmiginal
system and the identified model.

To my parents

ACKNOWLEDGEMENTS

I am very grateful to my thesis chairman and advisor, Professor Brian Feeny, for his
introduction to chaotic systems, for his understanding of my poor English, and his
encouragement and guidance in the research. His promptness and attention to detail were
extremely helpful. I thank Professor Alan Haddow for his excellent teaching in nonlinear
mechanical vibrations, for his friendship, and for helping my family in adjusting to a
foreign environment. I thank Professor Hassan Khalil, who also provided excellent
teaching in linear and nonlinear control theory, and stimulated ideas in my research. I also
thank professor Tien-Yien Li for participating in my thesis committee. My warm thanks
are also due to Professor David S.-Y. Yen, who left my thesis committee due to a medical
leave. I wish him the best for his health.

A special appreciation goes to Professor Joe Cusumano and Bart Kimble at Penn State
University for generously giving me the experimental data and detailed explanations of
the experiment. Without this, I would probably be working on this research for another
year.

A warm thank you goes to my colleague Shyh-Leh Chen, who was always enthusiastic
and eager to help in many ways. He provided many critical ideas and suggestions in my
thesis. Also, thanks to my other colleagues, Matt Brach, C. P. Chao, J. W. Liang, and
Ramana Kappagantu, for good discussions and companionship during my stay in MSU.

iv

Without my parents, I would not be here. My parents have never gone to school—not even
a single day in their lives-~but have worked very hard to send me and my brothers to
schoolforashighaneducadonuwecouldanaimmybelievethatacademic
achievement is the highlight of the family history. Their spiritual support played a crucial
roleinthecomplefionofthisthesialamalsoindebtedtomypuents-in-law, who
provided extra financial support for my family, subsiding my worries about living in
America. ’

My wife and my two daughters, Ting-fang and Yu-hua, came to this country shortly after
me. Ting-fang and Yu-hua came here with no English background. and were my main
oonoernotherthanmy study.'I'heyweottoSpartanVillageElementary SchooLwithgreat
courage and determination. They slowly but steadily improved their English, made friends
and enlarged their knowledge of American culture. Ting-fang went on to Hannah Middle
SchoolasaéthgradaandtheanacDonaldMiddleSchooluahhgradadmmthe
school district reorganization. Yu-hua continued through 4th grade in Spartan Village
School. They now enjoy their school lives very much. I am grateful to the faculty and Stafi‘
of Spartan Village School for pmViding mold-national cultural environments, and
excellent teaching for the international students. The same goes to Hannah Middle School
and MacDonald Middle School too.
OfcomsemywifeplayedanimportantroleregardinglifeinAmerica. Shetakescareof
ourdaily lives. She isavery goodcookthhout her, livinginAmericawouldheoome
more dimcult and less enjoyable. I deeply appreciated her moral support.

This research was supported financially by a grant from CSIST, Taiwan.

TABLE OF CONTENTS

1 Introduction 1

1.1
1.2
1.3
1.4

AnOverviewoftheParametricIdentificationMethods ............. 2
ChaoticMotion ............................................ 4
Motrvatron ................... 7
ThesisOverviews .......................................... 8

2 Methodology 10

2.1
2.2
2.3
2.4

2.5
2.6

Introduction .............................................. 10
Periodic Orbit Extraction .................................... 11
The Choice of a Mathematical Model .......................... l4
Algorithm ................................................ 17
2.4.1 Externally Excited Systems ............................ 17
2.4.2 Parametrically Excited Systems ......................... 21
2.4.3 Autonomous Systems ................................. 23
Strategy for Model Validation ................................ 24
Summary ................................................. 25

3 NumericalResults.........'............................. 26

3.1 The Forced Duffing Oscillator ................................ 27
3.1.1 Extraction of the Periodic Orbits ........................ 28
3.1.2 arousing a Mathematical Model ........................ 29
3.1.3 Identification Results ................................. 30
3.1.4 Model Verification ...... _ ............................. 31
3.1.5 Efl‘ect of Noise ...................................... 32
3.2 A Smooth Coulomb Friction System ........................... 3‘7

3.3

3.4

3.5

3.6

3.2.1 Effect of Noise ...................................... 41

A Parametrically Excited System ............................. 42
3.3.1 Effect of Noise ...................................... 46
An Autonomous System: the Lorenz Equations ................... 46
3.4.1 Effect of Noise ...................................... 51

A Case Study on Modeling the Nonlinearity with a Power Series . . . . 52
3.5.1 Using the Known Function in the Model .................. 54
3.5.2 Using the Power Series Approximation ................... 54
Conclusion ............................................... 56

4 ExperimentalResults................................... 58

4.1
4.2
4.3
4.4
4.5
4.6
4.7
4.8
4.9

Introduction ............................................... 58
Experimental Setup ......................................... 58
Phase-Space Reconstruction .................................. 60
Periodic-Orbit Extraction .................................... 64
Choosing a Mathematical Model .............................. 67
Data Processing Issues ...................................... 67
Identification Results and Model Verification .................... 69
Estimation of the Natural Frequency and the Damping Ratio ........ 74
Discussion ................................................ 76

5 ErrorsinParameterEstimates........................... 78

5.1 Introduction ............................................... 78
5.2 Errors Induced by Noise ..................................... 79
5.3 Errors Induced by the Periodic Orbit Extraction .................. 81
5.4 Sensitivity of the Parameter Estimates to Errors .................. 86
5.5 Using Several Periodic Orbits ................................. 88
6 ConclusionsandFutureWork............................ 90
6.1 Conclusions ............................................... 90
6.2 Future Work .............................................. 93

fi

BIBLIOGRAPHY ........................................ 95

LIST OF TABLES

Table 1 Identification results using individual periodic orbit for Duffing’s equation ....... 3 1

Table 2 Identification results using 4 periodic orbit data for Duffing’s equation .............. 32
Table 3 Identification results for Duffing’s equation using noisy data. ............................ 36

Table 4 Identification results using 4 periodic orbits for Coulomb friction system .......... 4 I
Table 5 Noise effect on identification results for Coulomb friction system ...................... 42

Table 6 Identification results using 4 periodic orbits for a parametrically excited system 45

Table 7 Noise effect on identification results for a parametrically excited system ............ 4 7
Table 8 Identification results for the Lorenz equation (a) ................................................. 51
Table 9 Identification results for Lorenz equation (b) ...................................................... 52
Table 10 Identification results for Lorenz equation with 1% noise ................................... 52
Table 11 Force and response in model (3.13) .................................................................... 53
Table 12 Identification results using the exact function in Model (3.15) ........................... 54
Table 13 Identification results using power series in Model (3.16) ................................... 55
Table 14 Identification results for the experimental system. ............................................. 69
Table 15 Comparison of the natural frequency and the damping ratio ............................... 7o

LIST OF FIGURES

Figure 2.1 (a) A close recurrence of a chaotic trajectory, and (b) a precise recurrence after

a gentle adjustment of the starting point of a chaotic trajectory ....... 12
Figure 3.1 Phase portrait of the Duffing oscillator ........................... 28
Figure 3.2 Some extracted periodic orbits of the Dufi‘rng oscillator .............. 29
Figure 3.3 The simulated chaotic attractor and some of the extracted periodic orbits of the

identified model ........................................... 33

Figure 3.4 Bifurcation diagrams of the Duffing’s equation (a) using the original equation,
and (b) using the identified model with the average values in Table 2 . 34

Figure 3.5 (a) A noise-free periodic orbit, (b) the noise-contaminated counterpart . . 35

Figure 3.6 The simulated chaotic attractor and some of the extracted periodic orbits of the

identified model using the noise-contaminated periodic orbits ....... 37
Figure 3.7 Bifurcation diagrams of the Duffing’s equation (a) using the original equation,
(b) using the identified model of the noisy case ................... 38
Figure 3.8 Phase portrait of a Coulomb friction system ....................... 39
Figure 3.9 Some extracted periodic orbits of a Coulomb friction system .......... 40
Figure 3.10 Phase portrait of the parametrically excited system .................. 43
Figure 3.11 Some extracted periodic orbits of the parametrically excited system . . . . 44
Figure 3.12 Phase portrait of a Lorenz system ............................... 48
Figure 3.13 Recurrence plot of the Lorenz system (period lengths are indicated in the num-
bers of time steps used in the numerical integration) ............... 49
Figure 3.14 Two extracted periodic orbits of a Lorenz system: (a) period length of 110
time steps, (h) period length of 144 time steps .................... 50
Figtne 3.15 Nonlinear function in a power series ............................. 56
Figure 4.1 Sketch of the experimental setup ................................ 59

Figure 4.2 Autocorrelation function of the experimental data ................... 61

X

Figure 4.3 Correlation function of the experimental data . . . . . l ................ 63

Figure 4.4 Singular system analysis of the experimental data ................... 64
Figure 4.5 Reconstructed phase space of the experimental system ............... 65
Figure 4.6 Reconstructed phase space in singular coordinates .................. 65
Figure 4.7 Some extracted periodic orbits from the reconstructed phase space ..... 66
Figure 4.8 Some extracted periodic orbits in singular coordinates ............... 66
Figure 4.9 Qualitative nonlinear function of the experimental system ............ 71
Figure 4.10 Bifurcation diagram of the identified model, (a) with the cubic nonlinear term
in the model, and (b) with the fourth power term in the model ....... 72
Figure 4.11 Phase portrait of the identified model with cubic nonlinearity ......... 73
Figure 4.12 Some periodic orbits extracted from the phase portrait of Figure 4.11 . . . 73
Figure 4.13 Transfer function of the right well of the experimental system ......... 75
Figure 5.1 Close look of the periodic orbit extraction on the Poincare section ...... 81
Figure 5.2 A sketch of the construction of a linearized map .................... 84

CHAPTER 1

Introduction

Parametric identification deals with the problem of determining the values of the
parameters in a mathematical model that represents a dynamical system, based on the
observed data from the system. It is a field of increasing interest, in part because of
applications in prediction and control.

Linear models have dominated in the description of the dynamical systems and control
theoretic approaches. Very complex and randomlike behavior has been viewed from a
statistical perspective in which many degrees of freedom were involved. More recently
nonlinear models have emerged, capable of describing chaotic dynamics and other
nonlinear behaviors. Such systems can exhibit extremely complex dynamical behavior,
even though the underlying dynamics may be low dimensional. On the other hand, high-
dimensional systems such as fluids and lasers can Show simple and low—dimensional
dynamics, which may be described by a low dimensional model [43].

Chaotic motion features the sensitive dependence on initial conditions. Nearby orbits that
cannot be distinguished will diverge exponentially and soon become uncorrelated. Along

with new theoretical concepts have come practical techniques, such as Lyapunov

2
exponents and fractal dimension, for characterizing the dynamics of such systems. Yet, the

techniques for identifying the parameters of a chaotic system are not as well developed as
those for analyzing its dynamical behaviors. We briefly review parametric identification
methods below.

1.1 An Overview of the Parametric Identification Methods

Parametric identification work generally presupposes that a mathematical model has been
chosen to represent a nonlinear system and that the goal is to identify the unknown
parameters in the given model. The unknown parameters are determined by optimizing in
some sense the fit of the chosen model to the available data.

For linear systems, the superposition principle can be applied to the system response and
the transfer function that characterizes the system behavior can be obtained by a variety of
techniques, such as transient analysis, frequency analysis, correlation analysis, and
spectral analysis [37]. The system parameters are estimated by a curve fitting of the
transfer function.

For a nonlinear system, the techniques for linear systems fail fundamentally because the
superposition principle is no longer applicable. However, for small nonlinearities,
perturbation techniques were widely used in analyzing the system response and in
identifying the system parameters as well. For example, Hanagud et al. [32] used the
method of multiple scales to determine the nominal system response, which was used
iteratively to estimate parameters. Nayfeh [45] and Feeny er al. [23] used the method of
multiple scales to exploit resonances and produce expressions relating the parameters to
the experimentally measured nonlinear behavior such as jump phenomena and nonlinear

drifi. The parameters could then be determined algebraically, or in a least-squares sense.

3
Ibanez [33] used a describing function method to construct an approximate transfer

function of the nonlinear structural system and hence uncouple the original nonlinear
equations. System parameters were obtained by iteratively minimizing the error function
between the measured data and the theoretical solution. Gottlieh et al. [27] used the
Hilbert Transform in their parametric identification of weakly nonlinear systems.

Another approach proposed by Mook er al. [41, 42, 52], called the method of minimum
model error (MME), combined the assumed model with the measurements to determine
the correct form of model for the nonlinear system under investigation. A correction term
which represented the model error was added to the assumed model and a cost function
was formed. By minimizing the cost function, a two—point-boundary-value problem was
formed and yielded the correction term, which was then fitted to an assumed polynomial
form to obtain the correct model of the nonlinear system.

Mohamrmd [40] used a direct approach by assuming a general form of the equation to
represent the nonlinear system under investigation. By measuring all of the system
responses, such as acceleration, velocity, and position, and directly introducing them into
the assumed equation, a set of algebraic equations was formed by balancing these
measurements and the input function. System parameters were then estimated by a
singular-value decomposition method.

In a similar direct approach, Yasuda et al. [59, 60, 61] represented the system nonlinearity
as a sum of polynomials in the system equation, with unknown coeficients as the system
parameters to be determined. Periodic responses under periodic excitation were measured
and expressed in Fourier series. The harmonic balance method was used to balance the

Fourier coemcients of each harmonic and a set of algebraic equations in the system

4

parameters was formed. The system parameters were then eStimated by a least squares fit
to the algebraic equations.

Recently, methods for modelling a chaotic system and identifying the system parameters
based on experimental data have been developed. Abarbanel er a1. [1] proposed a method
for constructing a parameterized map which evolved points in the phase space into the
future. This map was regarded as a dynamic system. and the parameter values were
chosen through a least-squares optimization procedure, constrained by the invariants of
the system, such as the Lyapunov exponents. Eisenhammer et al. [20] proposed a
trajectory method to extract ordinary difi'erential equations from an experimental time
series. The experimental data were represented in a state space and the corresponding flow
vectors were approximated by polynomials of the state vector components. Starting from
appropriately chosen initial states, the model equation was used to obtain an estimation of
the states for later times, and the coefficients were fitted by minimizing the distances
between the states predicted by the model and the experimental states. Breeden and
Hubler [6] proposed a flow method for reconstructing a set of coupled maps or ordinary
differential equations from a trajectory of the system in state space. By choosing some
trial coefficients for a series expansion in the state variable, the error in these parameters
were computed by comparing the predicted values and the experimental values. The

parameters of the model were obtained by solving a chi-squared minimization problem.

1.2 Chaotic Motion
Chaos was known by Henri Poincare (1854-1912) about a century ago in the course of his
investigations on the three-body problem. Through his discovery of homoclinic solutions

(homoclinic intersection, or homoclinic tangles), Poincare showed that the three-body

5
problem has no solutions of the type envisioned by Jacobi or Hamilton, in the sense that a

small error in the initial conditions produced an enormous error in the final response. In

his book, New Methods of Celestial Mechanics, Poincare motel:
When we try to represent the figure formed by these two curves and their
intersections in a finite number, each of which corresponds to a doubly asymptotic
solution, these intersections form a we of trellis, tissue, or grid with infinitely
serrated mesh. Neither of the two curves must ever cut across itself again, but it
must bend back upon itself in a very complex manner in order to cut across all of
the meshes in the grid an infinite number of times.
The complexity of this figure will be striking, and I shall not even try to draw it.
Nothing is more suitable for providing us with an idea of the complex nature of the
three-body problem, and of all the problems of dynamics in general, where there is
no uniform integral and where the Bohlin series are divergent.
In the words of modern dynamical systems dreary, the solution is sensitive to initial
conditions due to the inherent stretching and folding process of the nonlinear dynamics.
This sensitivity to initial conditions makes the nearby states on the attractor divergent
exponentially on the average, and results in a long-term unpredictability emanating from a
small amount of uncertainty in the initial conditions.
Confronted with his discovery of the homoclinic solution, Poincare went on inventing
several theories for new branches of mathematics, including topology, ergodic theory,

homology theory, and the qualitative theory of difi‘erential equations. He also pointed out

the possible uses of periodic orbits in characterizing his discoveryzz

 

1. “mummaanymuapawwmuwmommmcm 1992forrnore
historicalcornmatts[53].

6
there is a zero probability for the initial conditions of the motion to be precisely

those corresponding to a periodic solution. However, it can happen that they difi’er
very little from them, and this takes place precisely in the case where the old
methods are no longer applicable. We can then advantageously take the periodic
solution as first approximation, as intermediate orbit, to use Gylden’s language...
Given equations of the form defined in art. I 3 and any particular solution of these
equations, we can always find a periodic solution (whose period, it is true, is very
long), such that the difi'erence between the two solutions is as small as we wish,
during as long a time as we wish. In addition, these periodic solutions are so
valuable for us because they are, so to say, the only breach by which we may
attempt to enter an area heretofore deemed inaccessible.
The periodic orbits are dense in the chaotic attractor, and all of them are unstable. This is a
characteristic sign of chaos that only the presence of unstable periodic orbits but absence
of the stable ones [2]. The periodic orbit theme has been pursued by many authors in
modern dynamical system theory in characterizing a chaotic attractor [5, 19, 35, 49, 53],
and in the course of controlling a chaotic system [16, 47, 50]. The unstable periodic orbits
have also been used in system identification [31], and in recognizing parameter variations
[35]. We use them as a major tool in our parametric identification scheme for a chaotic
system.
Chaotic signals have been discarded in the past as ‘noise’. But, as pointed out by
Abarbanel [2], “chaos is not an aspect of physical systems which is to be located and

discarded, but is an attribute of physical behavior which is quite common and whose

 

2. MncKly, R. and Main, J.,eds., Hamiltonian dynamical systems (Adam Hilger, Philadelphia, 1987), cited from ‘l‘ufil-
m Abbott, an! Reilly inAnExperimental qpmaeh toNonlinearDynamicsandChaos.1992[53]

7
utilization for science and technology is just beginning’. It has been discovered in many

nonlinear systems in the laboratory and in the mathematical models in the past two
decades, and has become a well-known phenomenon and an important subject in modern
dynamical system theory. Much of the work has concentrated on learning how to classify
the nonlinear systems by analyzing the output from known systems. These efi‘orts have
provided, and continue to provide, significant insights into the kinds of behavior which
one might expect from nonlinear dynamical systems, and have led to an ability to evaluate
now familiar quantities such as fractal dimension, Lyapunov exponents, and other
invariants of the nonlinear systems [2]. Efl'orts have also been extended to predicting and
controlling the chaotic behaviors. For example, Farmer and Sidorowich [20] proposed a
local approximation approach for predicting a short-term chaotic time series using the
nearby states. Du et al. [47], Ditto et al. [16], and Shinbrot et al. [50] tried to control the
chaos by exploiting the periodic orbits embedded in a chaotic attractor and perturbing
some parameters of the system, so as to stabilize one of the unstable periodic orbits,
making the system become stable and more flexible under different operating conditions.
Cusumano and Sharkady [l2] experimentally studied the bifurcation and dimensionality
of a chaotic attractor occurred in a low dimensional parametric-excited system, and built a
valid model for the physical system. This trend of study shows that the chaotic motion

may be often regarded as an annoyance, yet it provides an extremely useful capability

without counterpart in non-chaotic systems.

1.3 Motivation
Chaos is inherent to nonlinear dynamical systems, and is rich in information content as

compared to a periodic trajectory. This richness has been exploited in dimensionality

8
studies, nonlinear prediction, and control schemes as stated in the previous section. The

potential usage of chaotic system in parametric identification has not been fully exploited,
because of the sensitive dependence on initial conditions and the long-term
unpredictability.

It is well known that a chaotic attractor is the closure of the set of unstable periodic orbits
[5]. They can be extracted and used to characterize the chaotic attractor [5, 19, 35, 49, 53],
and hence can be used for identifying the system parameters, because they are the solution
to the system equation.

Meanwhile, Yasuda and co-workers [59, 60, 61] have demonstrated that the stable
periodic solution to a nonlinear system can be used to identify the system parameters. This
inspires us to explore the applicability of the unstable periodic orbits in a parametric

identification scheme for a chaotic system.

1.4 Thesis Overview

In Chapter Two, we describe the methodology for identifying the parameters of a chaotic
system. A mathematical model is chosen to fit the characteristic of the original system
from which the chaotic data are obtained. The unstable periodic orbits are extracted from a
chaotic set for use in the identification algorithm. Then the harmonic-balance method is
applied to form a set of linear algebraic equations in system parameters, which are then
solved by a least-squares fit. This approach is applied to different kinds of nonlinear
systems, such as externally excited, parametrically excited, and autonomous systems.
Chapter Three contains the identification results for several numerical examples. Chaotic
data are generated numerically from known governing equations. Mathematical models

are chosen in polynomial form generally if no knowledge about the system nonlinearity is

9
available. Random noise is added to the periodic data to assess its efi‘ects on the

identification results. The identified models are verified by comparing the Lyapunov
exponents, the structure of the unstable periodic orbits, and the bifurcation diagrams of the
original system and the identified one.

In Chapter Four, we apply the method to a set of experimental data, taken from J. P.
Cusumano and B. W. Kimble at Pennsylvania State University. The phase space is
reconstructed by the method of delays, from which the unstable periodic orbits are
extracted for use in the identification procedtu'e. A model is obtained and verified for the
experimental system.

Two sources of error, noise and the extraction of the periodic orbits, in the identification
process are discussed in Chapter Five. We examine a bound on the error in the Fourier
coefficients induced by the noise. We examine the extraction of the unstable periodic
orbits closely, and establish a bound on the deviation of the extracted periodic orbit from
the real one. We discuss the sensitivity of the errors induced by the noise and the
extraction of the unstable periodic orbit to the identification results.

Chapter Six contains some conclusions and future work.

CHAPTER 2

Methodology

2.1 Introduction

Parametric identification method is not well developed for a chaotic system, partly
because in the past chaos has been treated as noise to be discarded, and partly because
chaotic motion exhibits sensitive dependence on initial conditions and defies long-term
predictability. Traditional usage of time series data in a parametric identification scheme
for non-chaotic systems may not be appropriate for chaotic systems, because of sensitivity
to initial conditions. However, a chaotic system features a chaotic attractorl, in which
infinitely many unstable periodic orbits are present, but absent of stable ones [2]. These
unstable periodic orbits can be extracted and used to characterize the chaofic attractor [5,
35]. They provide a skeleton of the chaotic set, which can be used in characterizing a
chaotic system.

Each unstable periodic orbit is a “solution” to the system which generated the chaofic set.
Once a periodic orbit is extracted, each term in the mathematical model can be expressed

in a Fourier series, and the harmonic-balance method can be applied to form a set of

 

1. Mmhmuficfingsetwhichconmadaueubinhisdiffiadtmshowmmpluthetaderrseorbit
Mmdmfactmmy oftlrenmnericallyobeavedfittrectorfmaymtbetrueattncm butmerely attractingsets.
Weusethistermlooselytodenoteasetofpoirrtsinplmespecetowudwhiehatirnehistoryapproechesaftertnn-
mdcmkGuckenhermenndHohnee [301andMoorr[43] forstrictmathernaticaldefinitionandexamphs.

10

11
algebraic equations, in which the system parameters are estimated by a least-squares

method.
The harmonic-balance method has also been used as a parametric identification technique
by Yasuda an coworkers [59, 60, 61] for some nonlinear systems that have stable periodic

response. We extend this technique to chaotic systems, where unstable periodic orbits are

abundant in the chaoric attractor.

In this chapter, we will demonstrate the methodology in detail with three kinds of chaotic
systems, categorized as externally excited, parametrically excited, and autonomous. They
are treated difl'erently because the excitation can affect the formulation of the
identification problem. We will discuss two important issues to our identification method.
the extraction of the unstable periodic orbits and the choice of a valid mathematical

model. We will also discuss the method for model verification.

2.2 Periodic Orbit Extraction

The genesis of a chaotic trajectory can be visualized as a random walk on the union of
infinitely many periodic orbits [14]. A physical trajectory approaches the saddle orbit
along its stable manifold, and remains nearby for a time before it is thrown out along the
unstable direction. It wanders around the union of periodic orbits, tracing out a strange1
attractor [14].

When the trajectory is near a periodic orbit, it approximately follows the motion of that
periodic orbit for an interval of time. If this time interval exceeds the period of the

reference orbit, the trajectory exhibits a recrurence. This property can be used to

 

1. mmnmmgemmfarmgmdBmmmmephnewmwhiehchmfieubhemehm
withageornetricobjectcalledefncteleet,whileadraoticefiractorfienotingebomdedmfionthntiemitiveb
changeehrhrifialeondifim,huuleastoneposidvel.yapmovexpmrau[43].Inouupecifictnrrpoee,weintendnot
todistinguishdwdifiamhnmdunmuchmgablymmfawdnm-tambehavbrofdnnmlmm

12
approximate the positions of the unstable periodic orbits embedded within the attractor [5,

35].

Figure 2.1 is a sketch of a recurrent three-dimensional flow in the vicinity of a hyperbolic
periodic orbit [53]. The chaotic trajectory has a recurrent segment, shown in Figure 2.1(a),
which is very close to an unstable periodic orbit, shown in Figure 2.1(b). We can gently
adjust the “starting point” of the trajectory segment so that the segment nearly coincides
with the unstable periodic orbit and returns almost precisely to its starting point [14, 53].
This idea has been used in a control scheme to stabilize one of the unstable periodic orbits

for a chaotic system [16, 47, 50].
In practice, we may have a sufficiently large chaotic data set {xi} , i = 1, ...N, in state

space. We scan the data set for recurrences by seeking points that come within a specified

spatial distance a of one another after a fixed elapsed time, such that [5, 35]

 

Figure 2.1: (a) A close recurrence of a chaotic trajectory, and (b) a precise recurrence
after a gentle adjustment of the starting point of a chaotic trajectory.

13
Putt-'19] $6, (2.1)

where K is the number of points per period of the unstable periodic orbit. If x‘. is a

recurrent point, then x‘. H, xi +2, , are likely to be “near” the unstable periodic orbit.

Thus, the segment of data, {x‘.,x‘.+ v ...,x‘.+K_1}, is then taken as the ‘approximated’

unstable periodic orbit. The ‘true’ one is generally unobtainable. The value i is the starting
point of the unstable periodic orbit, and is related to the phase angle associated with that
periodic orbit relative to the forcing function. It is important to record the phase angle for
later use in the calculation of the Fourier coefficients. This will become clear when we do
the calculation.

In a periodically forced system. all periodic orbits have a period that is an integer multiple

of the forcing period such that K = no, 2no, 3n0..., where n0 is the number of points per

forcing period [30]. However, in an autonomous system, such as in the Lorenz oscillator,
there is no such forcing period. Instead there are infinitely many unstable periodic orbits
with incommensurate periods. These incommensurate periods can be obtained using a
recurrence plot, which can be constructed by varying the period length and counting the
number of recurrent points found for each period length. The recurrence periods will be
clustered around certain values, indicating the periodicity of the periodic orbits and hence
the number of points in a period, which is then used as the fixed elapsed time K in Eq.
(2.1) for locating the periodic orbits, and also used as the fundamental period in the

calculation of their Fourier coefficients.

This procedure is quite successful in finding the unstable periodic orbits in many chaotic

14
systems. However, if the positive Lyapunov exponent associated with the unstable orbit is

large, then in one period the orbit will most likely have so departed from the unstable
periodic structure in phase space that one will probably not be able to identify the unstable
periodic orbit within the spatial criterion. In such case, the spatial distance criterion 8 can
then be relaxed to include more points in the neighborhood. For finding low-order
periodic orbits, say less than ten, it is adequate to set a to be 0.5% of the maximum extent
of the chaotic attractor [5, 35].

The searching process for the periodic orbits may reveal several distinct unstable periodic

orbits with the same period number. Nonetheless, all extracted periodic orbits can be used

in the parametric identification algorithm.

2.3 The Choice of a Mathematical Model

For the task of parametric identification, it is important to choose a valid mathematical
model to represent the physical system from which the measurements are taken. To do
this, we need to know the order of the system and the form of the system nonlinearity.
For typical mechanical vibratory systems, the system order is twice the number of degrees
of freedom. Also, for a nonlinear system to be deterministically chaotic, the system has to
be three or more dimensional. For a forced single-degree-of-freedom system, a general

mathematical model can be written as

mx+f(1.x.t) = 0. (22)
where the time variable is taken as an additional dimension. For an autonomous chaotic

system, a general mathematical model can be written as

Y = My) . (2.3)

15

forye R",n23.

In the forced case, the general function f (x, x, t) can take some specific forms when the
excitation and the nonlinearity of the system are known. For example, a system is
externally excited if it is modeled with an inhomogeneous term in the governing equation,
and parametrically excited if the system difi‘erential equations have time-varying
coefficients. Also, the form of the system nonlinearity can be determined using the
physical law that governs the system dynamics and the background knowledge about the
physical system. For example, sin(x) is usually used to model a pendulum system. A
power series can be used to model a system with an unknown smooth nonlinear function.
Whenever possible, models based on the physical mechanism should be employed. Thus,
for an externally excited single-degree-of-freedom nonlinear system, Eq. (2.2) can be

recast more specifically as
P
m! + 2 Bi,- (x. x) = E (t) . (2.4)
i I l
where E(t) is a known external excitation, fi (x, 2) are some known functions of x and x ,

p is the number of nonlinear terms in the system model, and m and B,- are the unknown

parameters to be determined. For a parametrically excited nonlinear system, Eq. (2.2) can

be recast as

r Pt
mx+ 2,2,0) { 2 Myrna} = o. (25)

8'81 jlll

where g ,- (t) are the known parametric excitation functions, fil- (x, x) are some known

16

functions of x and x , r is the number of excitation terms, p i are the numbers of the
nonlinear terms associated with each excitation functions, m and Bi]. are the unknown

parameters to be identified. And for an autonomous system, Eq. (2.3) can be written as

Pt

y,- = Z Bijhgj (”ti = 1, ..., n, (n 2 3) . (2.6)
fall

where y = [y], ..., yn] T, and hij (y) are the nonlinear functions of the state variable y,
and p,- are the number of terms in each equation of the model.

If the form of the system nonlinearity in Eq. (2.4) and (2.5) is unknown, but can be
assumed as a smooth function, and the system is operated in the neighborhood of the
equilibrium point, then the unknown function can be approximated by a truncated power
series. This is reasonable, because any smooth function can be represented by a power
series in some neighborhood of the origin (equilibrium). However, this approximation by
a power series may be accompanied by issues such as convergence and optimal truncation.
Ideas of convergence and divergence make sense when we consider infinite series. Since
we are using a truncated series, these ideas are not critical. If a power series indeed
converges to our function to be identified, then it is best to use as many terms as possible
without introducing numerical problems associated with large exponents. If the
underlying function has a divergent power series in the range of data, then there would be
some optimal truncation which is unknown. Thus, an imperfect identification result seems

to be the norm.

17
2.4 Algorithm

2.4.1 Externally Excited Systems

The mathematical model for an externally excited single-degree—of-freedom system is

chosen as Eq. (2.4)

P
mx+ 2‘, B,f,-(x.x) = Em. (2.7)

i I 1
where the external exciting function, E(t), is considered to be periodic with single

frequency to , such as

E (t) a aeos (art) . (2.8)
Upon extraction of the periodic orbits from the chaotic attractor, the nonlinear functions

become periodic and can be expressed in Fourier series, such as

xp (t) =a — 202;] {a cos(-k—‘) + b sing—i”) (2.9)
15,, (t) -.‘=. i (L?) {bjcos(‘Z-(£—t)—ajsin(J-—:”)} (2.10)
j- 1
it"E‘§.(’#2-°)’ijs(’—1”)+b8141—?)
fi(x,X)p2% 0‘r- 2 {cijcos owl?” til)+d sin(L(£—t)} (2.12)
j at

with the Fourier coefficients calculated as

 

 

 

 

 

 

mttl

balar,

(2.13)

(2.14)

dij=

cij =::%-.J: f(x,Jt) Pcos(— 1: )dt
2
7'

where the subscript p denotes the function being evaluated using the period-k data, T is the
period of the employed periodic orbit, and t) is the phase angle of the extracted periodic
orbit relative to the forcing. Since the phase angle has been included implicitly in the
periodic-orbit data and the nonlinear functions (the beginning of the periodic orbit is the
index of the phase angle), the limits of integration in Eq. (2.13) and (2.14) are used in the
numerical integration of the data. Ignoring the phase angle will cause an inconsistency in
the Fourier series representation in Eq. (2.9) to (2.12), and consequently produce incorrect
identification results.

Substituting these Fourier series into the model equation (2.7), and balancing the Fourier
coefficients of Eq. (2.13) through (2.14) of any set of harmonics, a set of linear algebraic
equations in system parameters can be constructed. This usage of the harmonic balance
method contrasts its usual usage for response analysis, where the ordinary difi‘erential
equation is known, and the efi‘ort is to solve a set of nonlinear equations in Fourier
coefficients. For systems forced with a single harmonic, and for autonomous systems, the
method of harmonic balance requires nonlinearity so that several harmonics can be

balanced.

In this thesis, we typically use the multiples of the primary harmonic. Thus, for the

wh

COI'

rcpn

matn‘

 

H2n

 

19

example of k = l , the balance equations, in matrix form, are

 

 

 

 

 

 

.- 0 COO Cm)1 ,. F 1
2 m.T 0
2 1
-0) b1 dll dpl B: ____ 0 (2.15)
2
_(nm) a. c1" c” LBP 3
2 -' l. ..
:(nm) bu d1" dP'U
01'.
At! = q, (2.16)

where a is the parameter vector of the system model; A is a (2n+1) x (p+l)

coemcient matrix, with each column containing the Fourier coefficients of the
corresponding term in the system model; q is a (p + 1) vector, containing the Fourier
coefficients of the external forcing function, which contains a non-zero element a in our
periodic excitation case; and n is the number of terms retained in the Fourier series
representation. For general values of k, the indices and frequencies in the elements of
matrix A are scaled by k.

If 2n = p and the matrix A is non-singular, the parameter vector u can be determined
uniquely. In practice, it is statistically better if the algebraic equation of Eq. (2.16) is
overdetermined, so that 2n > p . Consequently the exact solution will not generally exist,
but a best solution can be obtained by a method such as a least-squares fit. We seek a
solution that can minimize the average error in all of the equations. The error function is

most conveniently chosen as the sum of squares, or defined as [4]

20
e :2 IAa—ql. (2-17)
The solution that minimizes Eq. (2.17) is called the least-squares solution of the linear
system. The minimization of the error function is performed by setting the partial

derivatives of the squared error function with respect to the parameters equal to zero, i.e.

an/aot = 0 , which leads to the so-called normal equations as

ATMOI- q) = 0. (2.13)

and the least-squares solution of the parameters vector u is

—r
a = (ATA) ATq. (2.19)
Since the operation of a matrix inversion is less accruate and time consuming, the normal
equation is often not recommended in the numerical implementation. The most general

least-squares solution using the singular-value decomposition method is

a = Vz’rUT q = Aiq, (2.20)

where through the singular-value decomposition, A = UZVT , and AT is its pseudo-
inverse; U and V are the orthogonal matrices with each column consisting of the left and
the right singular vectors of matrix A respectively; and ET is the pseudo-inverse of 2 ,
which has the non-negative diagonal elements being the inverse of that of the
corresponding terms in 2. (See, for examples, Atkinson [4] and Strang [51] for a
geometric discussion).

Here arises a question as to how many terms should be retained in the Fourier series
representation of a periodic solution. Theoretically, the number of terms in the Fourier

series should be infinite, but Mickens [39] has shown that the upper bounds of the absolute

21
magnitudes of the harmonic coefficients decrease exponentially, such that they become

ineffective in the least-squares estimation procedures. We found that the rule of thumb for

retaining the number of terms in the Fourier series is

SSnSS. can
where n is the number of harmonics of the primary (driving) frequency. This limits the
number of unknown parameters in the model that can be estimated using a single periodic
orbit. However, we can use several periodic orbits to form several sets of algebraic
equations, thereby augmenting the matrix A to increase the redundancy of algebraic
equations for the least-squares estimation. This treatment can improve estimation result
even if the number of unknown parameters is not excessively large. Moreover, when the
parameter set is small, each set of algebraic equations formed by individual periodic orbit
can be used separately to obtain statistical information such as mean values and standard
deviations. This availability of several extracted periodic orbits fiom a chaotic set
increases the applicability beyond that of a simple periodic response, such as the case by

Yasuda and coworkers [59, 60, 61].

2.4.2 Parametrically Excited Systems

A parametrically excited system has time-varying coefficients in the governing equations
of motion. Examples 'of this kind of nonlinear system are a pendulum with a moving
support [46], a column with an axial time-varying force [46], and a flexible beam under an
electromagnetic force [12]. Previous studies have focused on dynamic stability, in which
the introduction of a small vibrational loading can stabilize (destabilize) a system which

was statically unstable (stable) [46]. Recent studies show that the system can exhibit

!f

char
The

cho

Sid:
{ac
div

aPP

Usi

COT.

Wit

Hen

thg.

22
chaotic behavior in a large range of parameters [12, 38].

The mathematical model for a single-degree-of-fi'eedom parametrically excited system is

chosen as Eq. (2.5),

r Pr
mx'+ 2 gi(t) {2 a... 0(a)} = o. (2.22)

i-l j-l

This model is a degenerate case of a parameter estimation problem because the right-hand
side of the equation is zero. Also, the Fourier series representations must account for the
fact that the system variables are coupled with a time-varying function. To proceed, we
divide through by m. The 1: term is taken as a known quantity by the fact that the
approximate periodic solution of the original system is known, and moved to the right-
hand side of the equation.

Using the extracted periodic orbits, the evaluated nonlinear functions in the model are
periodic. The excitation functions in time and the nonlinear functions in x and x are

combined together when they are to be expressed in Fourier series, such that

~ c.. "
8.30:.1. t) = g,o(t)f,-j (x. x) = 45+ 2 cijkcos(kmt) +dmsin (kart), (2.23)
t-r

with the Fourier coefficients calculated as

2 IT“-
2 +¢ (224)
dijt = 7‘ ¢ it,- (1,1,!) sin(lttnt)dt

Here, the phase angle is included in the combined nonlinear function g”. (x, x, I) through

the variable x. The limits of integration are chosen to match with the phase angle of the

CXt

Sui

COR

2.4

far:

We
Perle

Cant

pl’lasC

23
extracted periodic orbits.
Substituting the Fourier series into the model equation (2.22), and balancing the Fourier

coefficients of each harmonic, a set of algebraic equations such as Eq. (2.16) can be

constructed, and the parameters can be estimated by a least-squares fit.

2.4.3 Autonomous Systems

An autonomous system of dimension three or more can exhibit chaotic behavior [30]. The

famous example is the Lorenz equation, given by

Y1 = 0(Y2’yl)
Y2 = pyl‘yz-ylygg (225)
Y3 = "BY3+YIY2

There exists a ‘butterfiy shaped’ chaotic attractor, in some region of the parameter space,
which consists densely of infinite many unstable periodic orbits whose periods are

incommensurate.

To extract the unstable periodic orbits, the incommensurate periods have to be determined
by constructing a recrnrence plot, as stated in section 2.2.

The mathematical model for an autonomous system is chosen as Eq. (2.6). It can be

chosen more specifically if we know the type of the autonomous system under
investigation. In an experiment, each state variable y‘. must be measured. Using the
periodic orbits extracted from the chaotic attractor, each term in the model is periodic and

can be expressed in a Fourier series with the fundamental frequency as the one obtained

from the recurrence plot The Fourier coefficients are calculated as before, except the

phase angle can be ignored since there is no forcing function involved. Treating the y‘. -

the
me

anc‘

inv:
[11]
[27]
We

data
Data
The
tech,

knov

24

terms as known quantities, and balancing the Fourier coeflicients of each harmonic in each
equation, a set of algebraic equations of the form of Eq. (2.16) is formed. The system

parameters are then obtained using a least-squares fit.

2.5 Strategy for Model Validation

The final step, and perhaps the most difficult step, in parametric identification procedure is
the model validation. The objectives of validation are to seek answers to questions such
as: Is the identified model adequate? Under what conditions is the model representative
the system? Traditionally, the method of model validation is to simulate both system and
model under similar conditions and compare the respective responses. This is subjective
and lacks consistency for chaotic systems, due to the system’s sensitivity to the initial
conditions [3, 57]. More sophisticated criteria based on geometrical and statistical
invariants have been proposed, such as embedding trajectories [9, 43], Poincare sections
[11], bifurcation diagrams [3], Lyapunov exponents [1, 58], and the correlation dimension
[27], to characterize and compare reconstructed attractor and identified model.

We will seek consistency of the identification results from using different periodic-orbit
data sets. This is the most convenient way to check the quality of the identified
parameters.

The positive Lyapunov exponent is an invariant quantity of a chaotic system. Several
techniques have been developed into algorithm for estimating Lyapunov exponents from a
known dynamical system or from observable [5, 19, 35, 58]. We use the computer codes
by Wolf et al. [58] to calculate the Lyapunov exponents, which will be used in verifying
the identified model.

We will also compare the structure of the unstable periodic orbits that are extracted from

the l

The

chao

bylk

2.6

VVbli
unsui
cone:

physi

25
the original attractor and from the one generated from the identified model. This is

reasonable because the unstable periodic orbits are the skeleton of the chaotic attractor.
The structure of the periodic orbits would provide some geometric information of the
chaotic system that is useful to assess the quality of the identified model.

The criterion of bifurcation diagrams of the system and the identified model, as suggested

by Aguirre and Billings [3], will also be used as a supplementary criterion when available.

2.6 Summary

We have outlined a scheme for identifying the parameters of chaotic systems by using the
unstable periodic orbits that are extracted from the chaotic attractor. The method is simple
conceptually and easy to implement. Models are chosen based on the knowledge of the
physical system, or on approximation by a power series. Bach term in the mathematical
model is expressed in a Fourier series, and the Fourier coeficients of each harmonic are

balanced to form a set of algebraic equations in system parameters, which are estimated

by a least-squares method.

Methods for model verification are proposed.

This is an extension of an existing method, previously applied to systems with a stable
periodic response [59, 60, 61], to chaotic systems. By using the unstable periodic orbits,
the method exploits the structure of the chaotic set. Thus, we overcome issues such as

sensitivity to initial conditions.

Cl

N1

In th
math
to ex
cmbc
find:
set of

esu'm

ffiCtiu
talter:
Swen
“first
imm'a

53’s [cl-11

CHAPTER 3

Numerical Results

In the previous chapter, we presented an approach for identifying parameters in a
mathematical model of a nonlinear system that exhibits chaotic behavior. The strategy is
to exploit the chaotic attractor of the system and extract the unstable periodic orbits
embedded within it. The extracted periodic orbits are used to express each term in the
model in a Fourier series, and their coefficients of each harmonic are balanced to form a
set of algebraic equations in system parameters, which are then obtained by a least-square
estimation. This approach can be applied to a general class of nonlinear systems with
smooth nonlinearity.

In this chapter, numerical studies on the forced Duffing oscillator, a smoorh Coulomb
friction system, a nonlinear parametrically excited system, and a Lorenz equation, are
taken to demonstrate the applicability of this approach. Numerical integration of the
governing equations is carried out using a Stir-order Runge-Kutta method on a Sun
workstation. Typically, 50000 chaotic data points are generated for with a time step
interval of one-100th of the forcing period, or with 0.005 time step size in autonomous

systems.

3.1

The

It is

27
3.1 The Forced Duffing Oscillator

The forced Dumng oscillator is given as

mx' + cx + Bx + 1x3 =- acostnt. (3.1)
It is a classic differential equation that has been used to model the nonlinear dynamics of
mechanical and electrical systems. With B = 0, Eq. (3.1) is a model for a circuit with a
nonlinear inductor [55, 56], and with B<0,y>0, it is a model for the postbuckling

vibrations of an elastic column under compressive loads [44]. It can be written as a set of

first-order differential equations

11 = xz/m

xz = — cx2 - 311- 'yxl3 + acoscut (32)
to fit the format of the computer integration routine in public libraries such as IMSL.
This equation admits chaotic motions for a large range of parameters. We choose the
parametervaluesas m =1, c = 0.2, B = y =1,andtheforcingtermasa = 27,and
to = 1.33 [54]. These parameters are to be estimated by the present method.
Using the numerical data generated from the governing equation, a phase portrait is
constructed as shown in Figure 3.1. We see that the trajectory wanders around the phase
space in the attracting set. Any initial condition within the basin of attraction leads to the
same qualitative appearance in the phase space. This is the attracting set from which the
unstable periodic orbits are to be extracted.
Also, the Lyapunov exponents, indicating the average exponential rates of divergence or

convergence of nearby orbits in phase space, are calculated using the computer code by

Wolf at al. [58]. They converge to r, = 0.13, and x, = —0.468, indicating that the

28

 

15 7 Y T I T V T f

10

  
    
     

‘2‘ “If“

—10

 

 

 

‘P

j,
I
XOP
.e
n.-
al-
(I

Figure 3.1 Phase portrait of the Duffing oscillator

system of Eq. (3.2) is chaotic, since there is one positive value and the sum of them is
negative. (There is a zero exponent from the computer code, corresponding to the time

variable of the vector field. We omit it for convenience).

3.1.1 Extraction of the Periodic Orbits

The unstable periodic orbits can be extracted using the recurrence property of the chaotic
attractor, as stated in previous chapter. We repeat the idea here to emphasize its
importance. We scan the data set in state space forward to locate the recurrent points that

are close within a spatial distance a , such that

Pt + K—x‘l s c (3.3)
for a periodic orbit with K data points in the orbit. Here the index i is taken as the phase

angle of the periodic orbit relative to the forcing function, which is to be used in the

mmmn

woof

dx/dt

dx/dl

29
calculation of the Fourier coefficients. The spatial distance 8 is chosen to be 0.5% of the

span of the data set [5,35]. Some of the extracted periodic orbits are shown in Figure 3.2.

pd—1 pct—2

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Figrue 3.2 Some extracted periodic orbits of the Duffing oscillator

3.1.2 Choosing 3 Mathematical Model

To identify the system parameters, we need a mathematical model that can catch the
essential feature of the original system. Some a priori knowledge about the original
system will help choosing a valid model. In this case, we know that the system is an
externally excited, Duffing-type nonlinear system. Hence we choose a model in
polynomial form, which has been commonly used in modeling the Duffing type nonlinear

systems. The model with viscous damping is written as

30
p t'
mx'+a.x+ 2 Ba = acostut, (3.4)
i-l

where m, (x, and B5 are the parameters to be identified.

3.1.3 Identification Results

Applying the extracted periodic-orbit data to the mathematical model of Eq. (3.4), each
term in the model is expressed in a Fourier series. The phase angle associated with each
periodic orbit is included in the calculation of their Fourier coefficients, as discussed in
detail in Chapter Two. Then the principle of harmonic balance is applied to the primary
harmonics of the Fourier series, resulting in a set of algebraic equations in system
parameters to be estimated by a least-squares fit.

We first apply four sets of the periodic-orbit data separately to the model of Eq. (3.4), with
five terms retained in the polynomials, the identification results are shown Table l. ’
Also, we apply four sets of periodic-orbit data together to increase the redundancy of the
least-squares fit with different number of terms included in the model. The identification
results are shown in Table 2.

The results are accurate compared to the actual values, and consistent with each other for
using different set of periodic- orbit data. The non-zero parameter values are recovered
within 1% of their nominal values, and the zero-valued parameters are close to zero, even
when the mathematical model contains many unnecessary high-order nonlinear terms. The
standard deviations are less than 1% of the non-zero parameter values, or close to the
average values of the zero-valued parameters, indicating the consistency of the results.

Combining individual sets of algebraic equations increases the redundancy in the least-

 

square;
include

The id.

“0111111.

3.1.4
From t
the, id,
“MCI

“0313b

31

Table 1: Identification results using individual periodic orbit for Duffing’s equation

 

periodic

 

 

orbits (1 Bl 52 53 B4 as
actual 1.0000 .2000 1.0000 0.0000 1.0000 0.0000 0.0000
pd-3 0.9999 .1997 0.9915 -.0023 1.0009 .0002 .0000
pd-4 1.0008 .2002 1.0498 -.0050 0.9803 .0007 .0010
pd-5 0.9997 .1998 0.9659 -.0020 1.0124 .0001 -.0006
pd-6 1.0001 .1999 1.0015 .0016 1.0009 -.0002 -.0061
Average 1.0001 .1999 1.0022 -.0019 0.9968 .0002 -.0014
std. dev. 0.0005 .0002 0.0351 .0027 0.0134 .0004 .0032

 

squares fit, and improves the accuracy of the identification results when the model
includes many parameters.

The identified results suggest that the model can be refined by removing the higher-order
nonlinear terms whose parameter values are negligible. The reduction of the unnecessary

terms in the model tends to yield higher accuracy in the identification results.

3.1.4 Model Verification

From the numerical results, the model of Eq. (3.4) can be easily verified. However, we use
the identified model with the average values in Table 2 to generate a set of data, and
extract the unstable periodic orbits, for comparison with the original ones. The extracted

unstable periodic orbits from the identified model are shown in Figure 3.3. They resemble

32

Table 2: Identification results using 4 periodic orbitsfor Duffing’s equation

 

 

 

orders m (1 Bl 92 B3 94 95 55 97
actual 1.000 .21XX) 1.“)0 0.CXX) MIX) 0.000 0.000 0.1K” 0.(X)0
p84 1.1XX) 0200 1.005 -0.(X)l ”XXI 0.1K”

p=5 1.000 0200 1.015 p.001 0.996 0.000 0.000

p=6 MIX) 0200 1.010 0.(X)4 0.998 -o.oor 0.(XX) 0.000

ps7 l.(X)0 0.200 1.(X)5 0.004 1.000 -0.(X)l 0.(X)0 01X” -0.(X)0
Avg. 1.(X)0 0.200 1.009 .0015 0.999 -.(X)1 0.1K” 0.(XX) 0.(X)0
std.dv 0.000 0.000 0.005 0.(X)3 0.002 0.1K“ 0.(XX) 0.(XX) 0.1K”

 

closely their counterparts in Figure 3.1 and Figure 3.2.

The bifurcation diagrams, as shown in Figure 3.4, are calculated using the original
equation and the identified one, by slowly increasing the forcing amplitude and sampling
the steady-state response at the same Poincare section. The resemblance of the original
bifrncation diagram and the identified one can be clearly seen.

The Lyapunov exponents of the identified model calculated by the computer code of Wolf
et al. [58] are convergent to i, = 020, and 212 = —0.49, which are close to the original
values of 2.1 = 0.18 and 12 = -0.468 , with deviations of 11.1% and 4.7% respectively.

Thus our model is verified.

3.1.5 Effect of Noise

Numerically generated data are considered to be essentially noise free. The excellent

33

phone space period—1

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Figure 3.3 The simulated chaotic attractor and some of the extracted periodic orbits
of the identified model

identification results in the example above may deteriorate if noise is present in the
periodic data. To assess the influence of noise on the identification results, a set of
uniformly—distributed random noise is added to each periodic orbit for use in the
identification algorithm to test sensitivity of A, q and At: = q . If the noise is added to the
chaotic set before the extraction of the periodic orbits, the spatial criterion 8 may need an
adjustment.

The noise level is set by the ratio of its maximum amplitude to that of the employed
periodic-orbit data. Figure 3.5(a) shows a period—3 noise-free periodic orbit and Figure

3.5(b) shows its 2% noise contaminated counterpart in phase space. We examine the noise

 

 

 

 

 

 

 

 

 

 

 

as
tore. unplltudo

Figure 3.4 Bifurcation diagrams of the Duffing’s equation (a) using the original
equation, and (b) using the identified model with the average values in Table 2

eflect by using the model (3.4) with varying nonlinear terms and varying noise levels.
Four sets of noisy periodic-orbit data are used together in the identification algorithm. The
identification results are shown in Table 3.

Comparing with the previous results in Table 2 for the same model, we find that,

(1) within 5% noise level, the noise effect is not significant for a model with three or four
terms in the polynomial. As the nonlinear terms increase beyond five, the noise effect
increases. The errors are within 2.3% of the non-zero parameters. The zero-valued
parameters have larger deviations than in previous case;

(2) For the nonlinear terms in the model are greater than five, the noise efiect increase

rapidly, resulting in less accurate identification results. The issue of noise will be

35

 

 

 

 

 

 

 

 

 

 

 

,2
:0)
h
tel

Figure 3.5 (a) A noise-free periodic orbit, (b) the noise-contaminated counterpart

discussed in Chapter Five;

(3) For a small amount of noise, the effect of noise is not catastrophic to our method, but
in a robust way.

With the model identified using the noisy periodic data, we proceed to verify the model by
comparing the Lyapunov exponents, the structure of the unstable orbits, and the
bifurcation diagrams as before. Using a model with the parameter values as in the last

second row of Table 3, the Lyapunov exponents calculated by the computer code of Wolf
et al. [58] are convergent to £1 = 0.21 , and i2 = —0.5 , which are close to the original

values of )‘r = 0.18 and 3.2 = —0.468, with deviations of 16.67% and 6.84%

respectively.
The simulated chaotic attractor and the extracted periodic orbits from the identified model

Table 3: Identification results for Duffing’s equation using noisy data

36

 

 

noise orders m (1 Bl 52 B3 94 55
actual 1.000 .200 1.000 0.000 LNG 0.000 0.“

1% p=3 1.(X)0 .2“) 0.999 -.(X)6 MIX)
1% p=4 1.000 .2“) MIX) -.(X)1 1.000 .(XX)
1% p=5 1.000 .2“) 1.058 .016 0.983 .001 .001
2% p=3 1.1K” .2“) 0.999 -.013 1.(X12
2% p=4 1.000 .200 1104 -.029 1.001 .001
2% p=5 1.001 .201 0.947 .016 1.036 -.009 -.(X)5
3% p83 1.000 .2“) 0.995 -.020 1.(X)2
3% p=4 1.000 .2“) 1.002 -.042 1.1K” .001
3% p=5 1.(X)1 .201 1.136 -.044 0.958 .002 .002
5% p=3 MIX) .2“) 0.977 -.032 1.002
5% p=4 0.999 .200 0.988 -.069 1.001 .002
5% p=5 1.(X)1 .201 1.202 -.07 3 0.931 .(XB .(X13

 

37

 

 

 

 

 

 

 

 

 

 

 

 

 

“A an
—-v ‘v

--5 0 5 —5 0 I!
X X

 

 

Figure 3.6 The simulated chaotic attractor and some of the extracted periodic orbits
of the identified model using the noise—contaminated periodic orbits

are shown in Figure 3.6. The qualitative resemblance with the original ones of Figure 3.1
and Figure 3.2 is clearly seen.

A bifurcation diagram is constructed for the identified model, by slowly increasing the
force amplitude as the control parameter, and sampling the steady-state response at the

same time interval, as shown in Figure 3.7(b). It closely resembles the original one in

Figure 3.7(a). Thus the model is verified.

3.2 A Smooth Coulomb Friction System

A smooth Coulomb friction system is given as

x+cx+x+ (1+kx)tanh(ax) =fcos(tnt), (3.5)

“(111011 is one of the models of a dry-friction system, studied extensively by Feeny and

38

 

 

 

 

 

 

 

 

 

 

35
6 r r 1
E5 - . M
.. 0» 15K -
g3 / ..
as 210 l 3‘0 35

 

25
tomenplltudo

Figure 3.7 Bifurcation diagrams of the Duffing’s equation (it) using the original
equation, (b) using the identified model of the noisy case

Moon [24]. The dry-friction is often modeled as a multivalued, discontinuous nonlinear
force, which causes a “stick-slip” chaotic motion in a large parameter space. Here, the
smooth function, tanh (ax) , is used to approximate the Coulomb fiiction model. This
system exhibits ‘almost sticking’ motions, featuring a funnel-like structure in the phase
space under the harmonic excitation [24]. We choose the parameter values as c = 0.03 ,
k = 1.5, a = 50, and the forcing term as 1.9cos(1.3t) , for numerical simulation.
Numerical integration is carried out by a 5th—order Runge—Kutta method as before. A two-
dimensional phase portrait is shown in Figure 3.8, where a funnel-like structure is clearly

$6611.

The Lyapunov exponents are calculated from the known equation using the computer code

39

 

 
       
             

 

 

 

1.5 ~ -
1 t. ..
g 0.5 - , -
o - <r' ' ‘1 ... ’ ‘
.t¢//""-'3':' , ,
’35-'12»:
—o.5 - if," . " ~" -
-‘ I I A J;
-1 , .05 o 0.5 1 1.5

Figure 3.8 Phase portrait of a Coulomb friction system

by Wolf et al. [58], which converge to 711 = 0.13 and 12 = -39.8, indicating that the

system is chaotic indeed.

To proceed identifying the parameters of this chaotic system, we look for the periodic
orbits embedded within the chaotic set. Using the procedures stated previously, some of
the unstable periodic orbits are extracted, as shown in Figure 3.9.

Also we choose a model in a polynomial form as that of Eq. (3.4), assuming no knowledge
about the nonlinear function of the system. But, after conducting the identification

procedures, we found that the identification results were very poor. We postulate that the
power series representation of the nonlinear function, tanh (ax) , may not be valid with

the numerical data, due to the large value of a = 50. We then choose another model that

 

 

 

 

 

 

 

 

 

 

 

1 - 1.5 r
pd-l pd-2
0.5 » ‘ ’
0.5 . l
o t J
o P
—0.5 -0.5 i 1
'11) 5 0 05 1 :1) 5 0 0+5 1
2 f
pd-3
‘ i
o D
’11 0 1 2 '11 -o:5 0 0:5 1

 

 

 

 

 

Figure 3.9 Some extracted periodic orbits of a Coulomb fiiction system

contains the smooth function as a known function, such as

P . P .
mx’+ at + 2 aix' + 2 (Bix'-1)tanh (ax) = fcos (tut). (3.6)

is r i - 1
Applying four sets of the periodic orbits to this model with different numbers of nonlinear
terms retained, the identification results are shown in Table 4.
Up to the nonlinear order of five in the mathematical model, the parameters identified are
accurate within 1.0% error, with the standard deviation less than 1.0% of the non-zero
parameters or close to the average identified values of the zero-valued parameters.
The identified results suggest that the model can be refined by the same procedure as in

previous case. By removing the high-order nonlinear terms whose parameter values are

41

Table 4: Identification results using 4 periodic orbits for Coulomb friction system

 

 

 

P m c “1 0L2 as “4 Br 52 Ba [34
real 1.000 .0300 1.000 .(DOO .0(X)0 .(XDO 1.000 1.500 .0000 .0000
ps3 .9977 .0291 .9997 .0086 -.016 1.002 1.501 -.009

p=4 .997 1 .0294 .9904 .0192 .0036 -.021 1.002 1.503 -.006 -.008
p=58 .993 .030 .970 .030 .0170 -.021 1.000 1.515 -.000 -.036
avg. .996 .0295 .9867 .0193 .0015 -.021 1.(X)1 1.506 -.(X)5 —.022
std. .002 .0004 .0152 .0107 .0166 .000 .0011 .0076 .0046 .0198

 

a. a,=—0.003mda,=0.oo9

negligible, the accuracy of identification results is improved, and the consistency remains.
The model can be verified by the same criteria used in the previous case also. From the
numerical results that the identified values are close to the real ones within 1%, we expect

to get similar results, and reluctantly omit them here.

3.2.1 Effect of Noise

We wonder how noise will affect the identification results of this system. We add the
uniformly-distributed random noise to the periodic orbit for use in the identification
algorithm. Three sets of periodic-orbit data are applied to the chosen model of Eq. (3.6).
Parameter identification results are shown in Table 5.

We find that the noise deteriorates the accuracy of the parameter identification results

rapidly when the noise level is increased. Due to the large parameter value in the hyper-

42

Table 5: Noise effect on identification results for Coulomb friction system

 

 

noise level orders m c (11 0.2 0.3 B1 B2 B3
real 1.000 0.030 1.000 0.000 .(X)0 1.000 1.500 0.00
1% p=2 1.002 0.093 0.998 -.021 0.968 1.508
1% p= 0.994 0.165 0.96 -.005 0.933 1.515
2% p=3 0.996 0.084 1.019 0.017 -.083 0.979 1.518 -.041
2% p83 0.982 0. 150 1.004 0.072 -. 168 0.953 1.532 -.07 4

 

tangent function, tanh (ax) , the noise amplitude has been amplified significantly, such
that the accuracy of the identification results may have been distorted. This is a case in

which sensor noise may cause trouble in obtaining accurate identification results.

3.3 A Parametrically Excited System
The numerical example for a parametrically excited system is a nonlinear Matlrieu

equation, given as [13]

mx+cx+ (B+asint0t)x+u(y-asincut)x3 = 0. (3.7)
This is a model of an inverted pendulum under a two-well potential generated from above
by a magnetic dipole, studied experimentally by Cusumano and Sharkady [13] for the
bifurcation and low-order modeling of a parametrically excited system. The pendulum can

be buckled and unbuckled by changing the voltage applied to the electromagnet. It is also

an example of a periodically disappearing separatrix, analyzed previously by Coppola and

43

Rand [10], and Bridge and Rand [8].

The parameter values used to generate the chaotic data are m = l , c = 0.01 , B = 0.25 ,

‘y = 0.75, t) = 1/3, and the parametric forcing function is 0.55sin(0.28t) [13].

Numerical integration of the governing equations is carried out by a 5th-order Runge-

Kutta method on a Sun workstation as before. The phase portrait is constructed using the

numerical data, as shown in Figure 3.10, from which the unstable periodic orbits are

 

 

 

 

extracted, as shown in Figure 3.11.
0.: *fi
0.6- I.’ '3'», '2" fag-2.. '1 -1
e ‘.‘ ' :3 :3“.::" 0.
I: sx;§’,_zf-.r.f.ogf. ..
. .- equl'. “"-i§fi$ .- , xi
0,41- (‘.’:(r r) (1'; ,,:....’ g : ~1'.’2_ ..
e. (‘v .."p.‘r ‘% .u-A € i ' “it."
e ' a...“ . ‘ ' i ”s". 1' ~
5. .0 'V , I — ' ye \ . ‘5 :" .‘ \k ~
02' 3......{5‘ a! ‘21 r" . .rgfi““’i‘"}l’$‘i 7;:W-r-r-w'1's ‘
Susi... ' I J . . O '-:)‘:‘_;‘plfl(e.‘a'?;§5,&}‘ %, :‘g-Jz‘ fi )‘o‘
.‘ ... q.\ '3‘. .. g . ', .. .s .
01- 145W! : -3? r: Fifi-:3 "“551; “ we" ‘ ‘t‘ :‘r' -
_ “3.9).. 1 . 2‘. .91.“ {a ‘7' .‘315 """l§;'~ J
v' zoo ‘35!» 5."??? . «Tait-3° ' 3 art‘s-.3
.. :\.:.:c~..",: ,. _. J r g: Wfl‘é‘)’. ,‘ .r: “,‘ngr
.02. " '~.s- age - -6‘51- 3‘(~ ‘1""2r y, ”4, 1; -
I") o"m& . ' . ‘fl 1".“- i ' '1‘ ¢( ’9'}: 3
- ”if: '1 f ‘1' 1’ . .g, , ..g r- 3'. 3
4“ o-‘.‘:"'\?J.:;~’ tears as . arr- .» ...,.n. '
' ' ' ' r ’1, " r r ttr.-g.~=d‘r.;-: ‘-' . ' '
~ . .-:. ..-' J": ..rrsai"r' -‘ "
.05 . ~ ~ a," ”cigar-5e, -
.0.u I 1 A I I
—1.5 —1 —0.5 O 0.5 1
it

Figure 3.10 Phase portrait of the parametrically excited system

The Lyapunov exponents are calculated using the computer code by Wolf et al. [58], and

converge to A, = 2.03 and 2.2 = —32.35, indicating that the nonlinear system of

Eq.(3.7) is chaotic indeed.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

pd-1 pd-2
1 - 1 - -
0 5 r 0.5 r
g 0 r r g 0 r <
—O 5 1 .05 t
-1 ‘ ‘ * -1
-2 -1 0 1 2 -2 2
x
pd"!
1 - - 1
0.5 b 0.5
g o r g o l
v
-O.5 r -0.5 r <
-1 -1 ‘
-2 1 2 -2 -1 0 1 2
x

 

Figure 3.11 Some extracted periodic orbits of the parametrically excited system

The mathematical model is chosen as

P .
an 2 (Bi+yisin(ut)x' = —x . (3.8)

i=1

Here, the x-term is taken as an known quantity based on fact that the periodic solution is
known by the extracted periodic orbit data. Using these periodic orbits, each term in the
model is expressed in a Fourier series, with the parametrically excited force being
combined with the x-term as a single function. Applying the principle of harmonic balance
to each harmonic of the Fourier series, a set of algebraic equations in system parameters is
formed. Four sets of periodic-orbit data are used to increase the redundancy of the

algebraic equations, and the identification results are shown in Table 6, for difl‘erent

numbers of nonlinear terms retained in the model.

45

Table 6: Identification results using 4 periodic orbits for a parametrically excited
system

 

6 Bl $2 53 B4 55 71 72 73 74 75

 

 

real 0.01 0.25 0.00 0.25 0.00 0.00 0.55 0.00 - 0.00 0.00
.183

p=3 .011 .250 .001 .254 .550 -.00 -
.181

p=4 .011 .250 -.00 .252 .011 .550 -.00 - .010
.182

p=5 .011 .249 -.01 .249 .018 .012 .559 -.00 - .022 .034
.221

Avg. .011 .250 .004 .252 .015 .012 .553 -.00 - - .034
.195 .003

std. .000 .001 .006 .003 .004 .00 .004 .00 .021 .033 .000

 

The results are accurate compared to the real values for the number of nonlinear terms in
the model is not excessively large, such as the first two cases in the table (p=3 and p=4).
The non-zero parameters are close to the actual values, with less than 1% error of the
nominal parameter values, and the zero-value parameters are close to zero.

As the number of the nonlinear terms in the model increases, the identified parameters that
are coupled with forcing function are less accurate. Among the parameters, the damping is
small in value and vulnerable to numerical errors. In this case, perhaps the damping
should be estimated independently by traditional method, or using the global estimation
method proposed by Cusumano and Kimble [12].

Refining the model by removing the high-order nonlinear terms with negligible values, the

46
accuracy is improved. We expect the model verification will be satisfactory for this

parametrically excited system, as demonstrated by the example in section 3.1.

3.3.1 Effect of Noise

To assess the influence of noise on the identification results, the set of uniformly-
disuibuted random noise generated before is added to the periodic orbits for use in the
identification algorithm. The noise level is set as the percentage of the maximum
amplitude of the employed periodic data. Applying four sets of the contaminated periodic-
orbit data to the model of Eq. (3.8) with varying nonlinear orders, the identification results
are shown in Table 7.

Compared with the noise-free cases in Table 6, the efi'ect of noise on the identification
results are not significant for three or four nonlinear terms retained in the mathematical
model. It becomes more significant when high-order nonlinear terms are included in the
model and higher level of noise is added to the periodic orbits. Up to the noise level of 3%

of the employed periodic data, identification results are acceptable.

3.4 An Autonomous System: the Lorenz Equations

We take the Lorenz equation as an example of an autonomous chaotic system, written as

122 = px,--Jt2—.r11c3 (3.9)

33 = ‘ B’s +1112
with the parameter values as o = 16, p = 45.92, and B = —4. Numerical data are
generated from Eq. (3.9) using a 5th-order Runge-Kutta method for 10000 data points

with time interval being 0.025sec. There exists a ‘butterfly shaped’ chaotic attractor, as

47

Table 7 : Noise effect on identification results for a parametrically excited system

 

 

c B, B: B; B, BS 71 72 73‘ 7., 75
real .010 .250 .000 .250 .000 .000 550 .000 .183 .000 .000
1% p=3 .011 .250 -.00 .253 .547 .001 .177
1% p=4 .011 .251 -.01 .252 .01 .547 ..00 .177 .01
1% p=5 .011 .251 -.01 .243 .01 .02 .553 ..00 .208 .01 .03
2% p=3 .011 .251 -.00 .253 .544 -.00 .173
2% p=4 .011 .251 -.01 .252 .01 .544 -.00 .173 .01
2% p=5 .011 .253 -.01 .238 .01 .02 .548 -.00 .195 .01 .02
3% p=3 .011 .251 -.00 .253 541 -.00 .168
3% p=4 .011 .252 -.01 .252 .01 .541 -.00 .168 .01
3% p=5 .011 .254 -.01 .234 .00 .02 .543 -.00 .183 .00 .02
5% p=3 .011 .251 -.01 .254 .535 -.00 .158
5% p=4 .011 .252 -.01 .252 .01 .535 -.00 .158 .01
5% p==5 .011 .257 -.01 .227 -.00 .02 .533 -.00 .159 -.01 .01

 

a. negative vahrec as in the problem setting.

48

 

 

‘5
0

   

 

 

 

 

 

 

30,- 70- -
20- 60’
10- 50'

8‘! o- 940- ~
—10- 301' r
—20- 20> ‘
-301- 10~ r
4—340 -20 0 2‘0 40 340 -2.0 40 20 40

x1 x1

Figure 3.12 Phase portrait of a Lorenz system

shown in Figure 3.12, in which infinitely many unstable periodic orbits are dense and have
an incommensurate period associated with each periodic orbit.

Since the system is unforced, there is no fundamental period to be used as a guide for
finding the unstable periodic orbits for use in our parameu'ic identification scheme. We use
the recurrence property of the chaotic attractor to construct a recurrence plot to determine
the period length of the periodic orbits, as stated in detail in Chapter Two. The recurrence
plot is shown in Figure 3.13, in which the recurrent points that are clustered around certain
values can be clearly seen. These values indicate the incommensurate periods of the
periodic orbits. Using the values, the corresponding periodic orbits can be located within

the Lorenz attractor. Some of the extracted periodic orbits are shown in Figure 3.14, which

49

8

 

385'

§

T

144 .

220
259
* 330 352
168
158
101. 129
110
5- 56 l
.1 . . . . .
o 250

0 50 100 150 200 300 350 400 450 500
length of delay

 

 

momma“:
3 3 3

.5
0|
T

 

 

 

 

 

 

Figure 3.13 Recurrence plot of the Lorenz system (period lengths are indicated in the
numbers of time steps used in the numerical integration)

are to be used in our identification algorithm.
Knowing that the system is a Lorenz type autonomous system, a mathematical model is

chosen such that the linear terms and the quadratic nonlinear terms are included as

Jt1 = i (“1111+ ibif‘ixj] (3.10)

1.1 j21'

1-1 jzi

12 = i (“2131* £61}:ij (3.11)

3 3
x3 = (03.x, + Zdyxixj] , (3.12)
. 1

i jzi

50

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Figure 3.14 Two extracted periodic orbits of a Lorenz system: (a) period length
of 110 time steps, (h) period length of 144 time steps

where aij,

by, cij , and d". are the parameters to be determined.

Using the periodic orbits, each term in the model is expressed in a Fourier series with the
fundamental period obtained from the recurrence plot. The Fourier coefficients are
calculated as before, except the phase angle is ignored, due to the fact that there is no

forcing function involved. By balancing the Fourier coefficients of each harmonic in each

equation, and ueafin g the It i -terms as known quantities, a set of algebraic equations in
system parameters is constructed for the least-squares estimation.

We use two periodic orbits with period lengths of 110 time steps and 144 time steps
respectively, as shown in Figure 3.14, in the identification algorithm. The estimated

parameter values are shown in Table 8.

51

Table 8: Identification results for the Lorenz equation‘| (a)

 

1‘1 ‘2 J"3 11x1 xrxz 1113 J‘2“2 12x3 x3x3

x -15.951 15.971 -.098 -.030 .019 -.003 .000 .003 .003
1 (-16.0) (16.0)

45.748 -.927 .230 .039 -.028 -.995 .003 -.004 -.006
(45.92) (- 1.0) (- 1.0)

-.089 .041 -3.877 .031 .980 .003 .000 -0.003 -.004
(-4.0) (1.0)

a. Wvdmmtltcpar’arrtetervaluesofeachequationinfiq.(3.12)indicamdbythefiratcohnnn.

 

 

The actual parameters presented in the original system are identified accurately as
highlighted in bold-face in the table, although some of the zero-valued parameters are not
close to zero, such as the third term in the second equation.

The model equation of (3. 12) can be refined by knowing that there is no ‘square’ term in
the Lorenz equation. This refinement improves the accuracy of the identification results
significantly, not only the non-zero parameters are closer to the real values, but also the

zero-valued parameters are close to zero, as shown in Table 9.

3.4.1 Effect of Noise

To assess the influence of noise on the identification results, we add a set of uniformly-
distributed random noise to the periodic orbits as before. With 1% noise added to the
extracted periodic orbits, the identification results of the model without square-terms are
not significantly affected, although some of the zero-value terms have non-zero values, as

shown in Table 10. With higher-level noise added to the periodic orbits, the identification

52
results deteriorate rapidly. This case shows that noise is influential to the parametric

identification results for the autonomous system.

Table 9: identification results for Lorenz equation‘ (b)

 

Jr1 Jr2 Jr3 1132 11x3 Jr2"3

-16.011 15.953 -0.006 0.001 0.0002 0.001 1

 

x, (-l6.0) (16.0)

x 45.900 4.0228 -0.0127 0.0020 0.9996 0.0021
2 (45.92) (-1.0) (-1.0)

13 0.0041 0.007 -3.9991 0.9999 0.0001 0.0002

(-4.0) (1.0)

a. 111evalueaaretheparametervalueaofeachequationinfiq.(3.12)aaindieatedbythefiratcolunmwithoutthe
aquaretermainthemodel.

 

Table 10: Identification results for Lorenz equation' with 1% noise

 

 

 

x, J‘2 x3 x1x2 1113 Jr2"3
x -15.5813 15.9582 0.0757 0.0006 -0.0087 0.0011
1 (-16.0) (16.0)
1 47.4206 -1.8632 0.0256 0.0034 -1.0252 0.0112
2 (45.92) (-1.0) (-1.0)
x -0.9178 0.5717 -3.9468 0.9902 0.0200 -0.0129
3 (.40) (1.0)
a. Thevalueaaretheparametervaluesofeachequationinfiq. (3.12)asindicatedbythefirsteohum1withoutthe
square turns in the model.

3.5 A Case Study on Modeling the Nonlinearity with a Power Series

We have confronted a problem in modeling a hyperbolic-tangent function with a power

53
series in Coulomb—friction system. We postulate that the power series representation may

not be valid in such case. We examine this problem by a similar example, written as

1’ + cx + x + ktanh (x) = fcos (tut) (3.13)
We want to show that, if the nonlinear function is known, our method is capable of
identifying the parameters accurately, as shown in section 3.2; if the nonlinear function is
unknown, and the power series is used to approximate it, then the radius of convergence of
the power series and the truncated series representation are the factors influential to the
identification results.
The parameter values in Eq. (3.13) are chosen as c = 0.3, k = 0.5, and to = 1.3.
Numerical data are generated using the Runge-Kutta method with several forcing
amplitudes. The maximum periodic responses under difi'erent forcing amplitudes are
listed in Table 11.

Table 11: Force and response in model (3.13)

 

 

 

 

 

 

 

case force, f max. x
a 0.1 0.25
b 0.5 1.0
c 1.0 2.0
d 2.0 3.0

 

 

Note that, by Taylor series expansion, tanh(x) can be represented by

_ _13 2 s_1_7 7 3
tanh(x) -x 3x +15x 3151c +...,Ix1$2. (3.14)

54
3.5.1 Using the Known Function in the Model

We choose a mathematical model in a polynomial form as that of Eq. (3.4), with the

addition of the known function of tanh(x) in the model, such that

p .
mx+ax+ 2 Bix'rytanh (x) = fcos (cot),

i-l

(3.15)

where the parameters m, a, B: and ‘y are to be determined using the periodic data. The

identification results are very accurate, as shown in Table 12, even when the mathematical

model includes many unnecessary terms.

Table 12: Identification results using the exact function in Model (3.15)

 

 

 

cases 77: a 151 152 53 B4 ' BS 7

actual 1.0000 0.3000 1.0000 0.0000 0.0000 0.0000 0.0000 0.5000
a~b 1.0006 0.3000 1.0300 -.0009 -.0000 0.0000 0.0021 0.4700
a~c 0.9999 0.3000 1.0008 -.0000 -.0002 0.0000 0.0000 0.5000
a~d 1.0000 0.3000 1.0000 0000 .0000 0.0000 0.0000 0.5001

 

3.5.2 Using the Power Series Approximation

Assuming that the nonlinear function of the system is unknown, our first choice is using a
power series to approximate it. Part of the reason is that it is “easier” to fit the nonlinear
function with a polynomial, and “possible” when data is within the radius of convergence.

A model is chosen in a polynomial form as

55

P .
mx+wz+ 2 Bax‘ = fcos(tot). (3.16)
i=1

Applying the periodic data to this model, the identification results are liable to errors,
depending on the amplitude of the response and the nonlinear terms retained in the model,

as shown in Table 13.

Table 13: Identification results using power series‘ in Model (3.16)

cases 771 a 151 B; 133 B4 35 136 B7 138 B9

a~b 1.01 .300 1.52 .000 -.155 -.000 .036

a~d .951 .300 1.35 .000 -.062 -.000 .003

a~b .999 .300 1.50 -.000 -.166 .000 .060 -.000 -.014

a~d .991 .300 1.46 -.000 -.103 .000 .013 -.000 -.001

a~b 1.00 .300 1.50 .000 -. 167 -.000 .067 .000 -.025 -.000 005

a-d 1.11 .300 1.51 .000 -.139 -.000 .031 .000 -.003 -.000 .000
000

actual 1.00 .300 1.50 .000 -.l67 .000 .067 .000 -.027
a. Themnsretamedhrtheserieaiahtdicatedbythelstcolmnnmmber.

 

 

 

 

 

 

In each case, better identification results are obtained using the smaller response data
(cases a and b), which are within the radius of convergence of the series. The best result is
obtained in the last case, in which the smaller response data are used in the model with
nine terms included, which almost fits the power series in Eq. (3.14).

Although the nonlinear function of the system is unknown, we may obtain a qualitative
feature of the nonlinear function from the identification results. The nonlinear function is

plotted using the identified values, as shown in Figure 3.15, in which the qualitative

56

0.5 I I V I r T T

 

0.4 - ad-5 is the curve using the values in the case of a~d with p=5.

  

0,3 . Analogous to other curves.

 

  

 

 

 

0.2 -
0.1 - ' if i
3 o - .
-o.1 «
-o.2 .1
-o.3 -
-0.4 ab-S, a“ ab—7, ab—9 , 0.5mm -
"0'31 —o:a -o:e —o:4 —o:2 1:: 0:2 0:4 0:3 ate 1

Figure 3.15 Nonlinear function in a power series

feature of the nonlinear function is clearly seen.

3.6 Conclusion

Numerical examples taken from the Duffing’s equation, a Coulomb friction system, the
nonlinear Mathieu equation, and the Lorenz equation, show that the present method can
accurately identify the parameters in a mathematical model that has been well-chosen to
match the characteristic of the original chaotic system.

The mathematical model can be refined by removing the unnecessary terms that have
negligible values. Consistent identification results are remained for the valid models,
implying that the suspicious terms are indeed unnecessary. Models are verified by
comparing the structure of the unstable periodic orbits, the Lyapunov exponents, and the

bifurcation diagram. The usage of many periodic orbits in the identification scheme

57
improves the accuracy of the least-squares estimation and provides the statistical

information of the identified results. This is suitable for systems with many parameters to
be identified.

Random noise added to the periodic orbits can deteriorate the accuracy of the
identification results, but in a robust way. This efi'ect worsens when the mathematical
models are not well-chosen, for example with many unnecessary terms.

When the precise form of the nonlinearity is unknown, yet smooth, the accuracy of
identified truncated power series coefiicients deteriorates. However, the truncated power

series may be applicable for qualitative modeling.

 

 

CHAPTER 4

Experimental Results

4.1 Introduction

In this chapter, we investigate a chaotic data set taken from a periodically driven magneto-
oscillator by J. P. Cusumano and B. W. Kimble at Pennsylvania State University. The
experiment was designed for observing the global phase—space structure of basins of
attraction and homoclinic bifurcation using the stochastic interrogation method [12]. The
experimental system was known to be similar to a two-well potential system.

The techniques developed in the previous chapters are to be applied to the given set of
chaotic data, in effort to identify the parameters of this experimental system. The chaotic
attractor is reconstructed using the method of delays [26, 53], from which the unstable
periodic orbits are extracted for use in the identification algorithm. A mathematical model
is chosen in polynomial form by knowing that the experimental system has a smooth two-
well stiffness potential. The method of harmonic balance is used to form a set of algebraic

equations in system parameters, which are estimated by a least-squares fit.

4.2 Experimental Setup

The experiment conducted by Cusumano and Kimble consisted of a stiffened beam

58

59
buckled by two magnets. The beam had extra rigidity in the form of steel bars epoxied and

bolted along the length away from the clamped end. This additional rigidity was included
in effort to make the system behave as a single degree of freedom. The uncovered portion
of the beam near the clamped end acted as an elastic hinge from which the position of the
beam was measured by a strain gauge. Two rare-earth permanent magnets were placed on
the base of the frame holding the beam to create the two—well potential. The frame was
then fixed through a rigid mount to an electromagnetic shaker. A periodic driving signal
was fed through a power amplifier to the shaker to provide the external forcing function.

The experimental set-up is shown in Figure 4.1. [12]

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

A .1.
_ «MM 0-‘9 sauna-n.
'1 T‘—
I: 3 Ho
E 5 7-“ JOF
E Mao- O
3/ 3,. _l
2.11% J

 

 

 

.1 I... |._ *|°-’|'—

Figure 4.1 Sketch of the experimental setup.

Data from the strain gauge was acquired using a 12-bit data-acquisition (AID) board, with
the digital values from -2408 to 2407 corresponding to -5V to 5V. With no forcing, three

equilibria exist; two are stable at digital values of -495 and 315, and one (saddle) is

60
unstable at approximately zero. When forcing is added, periodic orbits exist instead of

equilibrium points. The driving frequency was set at 7.5 Hz with 1.5V of the function
generator output, by which the chaotic data were generated and collected at the sampling

frequency of 187.5 Hz for 7000 periods of excitation.

4.3 Phase-Space Reconstruction

Since there is only one observable in the data set, denoted by {xi} , j = 1, ...N, with

1]. = xUAr) , A: is the sampling time interval, the phase space of the experimental

system is to be reconstructed. The most common method of phase space reconstruction is
the method of delays [26, 53]. It is used to construct a d—dimensional pseudo-vector with

its elements being the single observable separated by a constant delay time, such that

yj = (xi’xjn’ ""xj+t(d-1))’ (4.1)
where 1: is the delay time, and d is embedding dimension. Both of which are to be
determined. The pseudo-vector represents a data point in the embedding space.

In theory, for any sufficiently large dimension d and almost any choice of delay time t , an
embedding of the original attractor can be obtained, and the geometrical invariants such as
dimension and positive Lyapunov exponents can be preserved. In practice, the delay time

I should be chosen so that the elements of yj are uncorrelated. If t is too small, then the
coordinates at X]. and xj +1 represent almost the same information. If t is too large, then

xj and xi +1 represent distinctly unrelated components of the embedding space. If the

embedding dimension dis too small, the trajectory may cross itself. The requirement of a

61
sufficiently large embedding dimension prevents such ambiguity and ensures that the

reconstruction is difl‘erentiable and invertible [2]. But an excessively large embedding
dimension may lead to excessive computation and corrupt data, since noise will dominate
the additional dimensions of the embedding space where no dynamics is operating [26].

There are several methods that have been proposed to determine the suitable delay time
and the embedding dimension [1, 5, 9, 25]. We use the criterion proposed by Abarbanel
[1] to determine the delay time 1: to be approximately 1/10~ 1/20 of the time
associated with the first local minimum of the autocorrelation function of the

measurement data {xi}. The autocorrelation function is defined as

1 N
R“) = [(7211431) (43)

ill

and is shown in Figure 4.2 for the chaotic two-well data.

 

 

 

 

 

O 20 4o 60 80 100 120 140 160 180 200
time stop.

Figtu'e 4.2 Autocorrelation function of the experimental data.

The pro;
the sing

015011»

whcr

then
in t‘
[2.

(1%)

cl:

62
The proper embedding dimension is estimated by the correlation function method [27] and

the singular system analysis method [9]. The correlation function method calculates the

distribution of points within a small region for a large data set, such that

N N
C(r) = P30 5. 2 2410-14-81) <44)
3 i

where H(z)=l if z is positive; and H(z)=0 otherwise; y is a pseudo-vector constructed as
Eq. (4.1). If the attractor is properly unfolded by choosing a sufficiently large dimension,
then any property associated with the attractor which depends on distances between points
in the phase space would become independent of the value of the embedding dimension

[2, 15, 27]. In a regime that C(r) becomes independent of d, and exhibits a power law

dependence on r as r -> 0 , that is limoC (r) = ard, the correlation dimension could be
r —9

obtained by measuring the slope of the plot of logC(r) versus logr, such as

d = lim M (4.4)

r -4 0 logr
Figure 4.3 shows the plot of logC(r) versus log(r) for several values of the dimension d.
The slopes are about 2.5, which becomes independent of the dimension as d 2 3 . Ding er
al. [15] reported that the plateau begins when the embedding dimension first exceeds the
correlation dimension. Thus, this criterion should produce a lower bound to our required
embedding dimension.
The singular-system analysis method involves constructing a covariance matrix C = YT Y

and decomposing it into two unitary matrices U and V and a diagonal matrix 2 , such that

$08.

Whe

Sin E

She

63

Correlation function mothod

 

 

 

 

Figure 4.3 Correlation function of the experimental data.

C = UZVT , (4.5)
where Y is the mauix with each column containing the pseudo-vector y as constructed in
Eq. (4.1). Varying the dimension in constructing the pseudo-vectors, and conducting the
singular values analysis, a plot of the singular values versus the embedding dimension is
shown in Figure 4.4.

By comparing the singular values with the values induced by ‘noise’, which is assumed
uniformly distributed in the extra dimensions and will be nearly equal, the singular values
become flat when d 2 4 . Thus we determine the suitable embedding dimensions to be
four.

A two dimensional projection of the reconstructed phase space is shown in Figure 4.5,

from which the unstable periodic orbits are to be extracted for use in our identification

-. 829 3:9...

algo

4.4

Fro

des

W1".

Ex

64

Singular value. vorouo ombodding dimension.

 

 

 

 

éau: d=6 12
400.. ............... § .............................. § ................ , ............................... _
35° .............................. Hd= .................................... -
30°. ...... oood=4 ................................... .1
325°). ........................................... §.m..®3 ................................... .1
3,..- ........................................... ................................................ _
“,0, ............................................ ............................................... .
m. .......................................... ................................................ ..
,0. .......................................... ............................................... -
oo 1 2L 4:_ 5 6
dimonoiono

Figure 4.4 Singular system analysis of the experimental data.

algorithm. Also, a uansformation of the reconstructed phase space into the singular

coordinates is performed by using the singular vectors, as shown in Figtue 4.6.

4.4 Periodic-Orbit Extraction
From the reconstructed chaotic attractor, the unstable periodic orbits can be extracted as

described in Chapter two. In the pseudo phase space, we seek recurrent points such that
lynr’yi $8 (4.6)

where e is set as 0.5% of the maximum extent of the chaotic set as before. Some of the

extracted periodic orbits are shown in Figure 4.7. The corresponding periodic orbits in the

singular coordinate are shown in Figure 4.8.

65

 

Reconetructed pheee epece of the experimental eyetem

2 f T Y Y

 

 

 

 

 

 

 

 

 

 

1 .5 - ~
1 r- ..
0.5 - ~
4 .— .
S’
—0.5 r ..
_1 I- d
-1 .5 " -4
—2 —1 5 —1 -0.5 0 0 5 1 1 5 2
X“)
Figure 4.5 Reconsu-ucted phase space of the experimental system
Projection of the reconetructed pheee epeoe on einguier ooordinetee
‘ .5 U V I I 1 V
1 .
.,v ’ ;~.’.-~-‘...\ r.‘\~ .
g 0.5 _ / I. _ — x
‘ " ,”l ’5 7//, .‘..
/ [7({1'T4g‘pti W\\ \\ ‘\
.. ' «((10.01
. ' \ \ t “ n 1/ .‘
g ‘ ‘- \\.:\‘\ \\’;.‘:.?.’/,4/,/ {/11 '
\ ::\\<f/ " '
\ .\~_ - g . ,
.0.‘ " \\‘\‘ t“::';:; ,.
-1 r-
-‘ .5 L 1 l l l l
—4 —3 -2 —1 0 1 2 3

eingular coordinate. q1

Figure 4.6 Reconstructed phase space in singular coordinates

38...: 33::

pd—1 psi-2

 

 

11W)
9

x(i+hu)
9

 

 

 

 

 

 

 

 

'32 o i 2 '52 -i 6 1 2
11(1) Xi!)
psi-43 psi—4
2 - - - 2 5 - -
j p 4 1 1

41M
9

#14181)
9

 

 

 

 

 

O 1 2 32 -1 0
x10 xtt)

 

 

Figure 4.7 Some exuacted periodic orbits from the reconstructed phase space

pd-1 pd—2

 

 

 

 

 

 

 

 

-1 .5 7 7 7 -1 .5 7 7 7
-4 -2 O 2 —4 —2 0 2
einouiar coordinate. q1 einguiar coordinate, q1
pd-s rad—4

 

fi v v

3?
0.5» « 0.51 , .
., O ‘ 3 °~ 2 ‘
-O.5r . —o.5r ‘ .
i» « i-u .

k ‘

 

 

 

 

 

 

 

-1. 7 - . 7 r 7
E4 —2 0 2 ‘ E4 —2 0 2
einguiar coordinate, q1 einguiar coordinate, q1

 

 

Figure 4.8 Some extracted periodic orbits in singular coordinates

67
4.5 Choosing a Mathematical Model

Knowing the experimental system is an externally excited nonlinear system with a two-
well potential, we choose a mathematical model in a polynomial form to fit the
characteristics of the nonlinear function. We choose a polynomial because we know that
the magnetic and elastic forces are smooth. We do not know, however, whether a power
series converges to the actual stifi'ness characteristic in the domain of displacement.
Furthermore, in the case of divergence, we do not know the optimal truncation of the
series representation. Our best hope is to obtain a model which qualitatively fits the
characteristic of the experimental system.

The model with viscous damping is then written as

P .
mx + at + 2 [53' = acosror, (4.7)

i- 1
where 77:, 0t, and B, are the parameters to be determined, p is the number of terms in the

power series.

4.6 Data Processing Issues

The experimental data are in a digital format, ranging from -2048 to 2047 corresponding
to -5V to 5V of the voltage output from the AID converter. There is a scaling factor
between the digital numbers and the actual physical unit. The parameters in Eq. (4.7) are
scaled by this factor in a nonlinear fashion. Assume that the factor between the digital data
2 and the variable x in Eq. (4.7) is a constant 7 in units of (displacement unit)/(digital

unit), such that

68

x = 72. (4.8)

Substituting this into Eq. (4.7), the model equation can be rewritten as

P . .
(my) 2+ (cry) 2 + 2 (1357‘); = acostut. (4.9)

i- 1
Since the digital data are large in amplitude, the high-order nonlinear terms will be even
larger in amplitude, causing an ill-conditioning of the matrix A used in the least-squares
fit. To prevent this, we can choose 7 in such a way as to normalize the data to the unit

interval. The time variable can also be nondimensionalized to a new variable, 1 = (or.
This normalization of time is manifested in the velocity and acceleration terms, and
improves the conditioning of the least-squares problem.

Meanwhile, we know that the external forcing is periodic, although the forcing amplitude
is unknown. This implies that Eq. (4.9) is actually indeterminate, and one of the quantities
in the equation has to be taken as known and moved to the right hand side of the equation,
as having been done in the parametrically excited and the autonomous cases in previous
chapter.

In this work, we deal with the unknown forcing amplitude by discarding the first harmonic
of Eq. (4.9) so that the other parameters may be identified. Another approach would be to

include the forcing amplitude as unknown sine and cosine coefficients to be identified.

Hence, we divide through Eq. (4.9) by the quantity myth2 , and recast it in a fonu as

1» ~ .
612+ 2 [3,2' = z"+2icosr, (4.10)
1-1

69
where (i = a/mtu, 8,- = Biyi-l/mtuz, and ii = a/m‘ymz, arethe scaled parameters to
be determined.
Using the extracted periodic orbits, each term in the model equation (4.10) is periodic and
expressed in a Fourier series as done before. The Fourier coefficients of the multiples of
primary harmonics, except the first harmonic, are balanced to form a set of algebraic
equations in system parameters for least-squares estimations as usual. The phase angles

associated with the extracted periodic orbits are ignored since the harmonics which

balance the forcing function are not used.

4.7 Identification Results and Model Verification
Using ten extracted periodic orbits together in the identification algorithm, with the data
being processed as discussed above, and using four terms in the polynomial in the model

Eq. (4.9), i.e., p=4, identification results are shown in Table 14.

Table 14: Identification results for the experimental system

 

 

 

p31,? a 6', 6'2 63 ii.

1 to 10 0.034 -0.266 0.141 0.323 -0.041
6, 8, 10 0.025 -0.304 0.184 0.338 -0.07 8
3, 5, 7, 9 0.041 -0.269 0.190 0.346 -0.067

3 to 9 0.034 —0.244 0.199 0.338 -0.070
Average 0.035 -0.280 0.178 0.336 -0.064
Std. dev 0.007 0.041 0.026 0.0096 0.016

 

The idcr
orbits. '
0011351
factor
We h:
stiffm
mega
dam
eno

As

let
Tc
fu

a\

[E

h4_2
—

70

The identification results are consistent by using difi'erent combinations of the periodic
orbits. The standard deviation of each identified parameter is small compared to the
corresponding average value. The values in the table are scaled to SI units by an unknown
factor as discussed in previous section.

We have also used a third-degree polynomial in the identification process. The identified
stiffness parameters were similar to those in Table 14. However, the damping term was
negative. A priori knowledge tells us that the damping should be positive. Since the
damping is small, it is likely that slight inadequacy in the third-degree model caused
enough error in the damping term to reverse its sign.

As in Chapter three, we are tempted to neglect the 4th-degree term since its coefficient is
small. However, we favor the identified parameters based on the 4th-degree model since it
led to a reasonable damping term.

To show effect of the 4th-degree term in range of data, we first examine the nonlinear
function of the system by plotting it with the identified parameter values. Using the
average values in Table 14, the nonlinear function is shown in Figure 4.9.

Within the scaled data range, the curves are qualitatively similar using three and four
terms in the power series. The curves represent a nonlinear function similar to the one
generated by a two—well potential to which the experimental system belongs. Hence, we
obtain a qualitative model for the experimental system, with unknown factors as discussed

in previous section, Eq. (4.10), in the following equation:

t' + 0.0352' — 0.282 + 0.17322 + 0.33623 — 0.06424 = am}. (4.11)

We proceed to do a numerical simulation using the identified model of Eq. (4.11).

However, since the force amplitude is unknown, we estimate it indirectly. By substituting

71

 

 

 

 

‘3 Y T Y T I
1 - ............... p:ua:..-o.za°z}o.1zarzezsogaaa'zea-rooeyrzu .......... .............. ii
03 t ............... fem ;.aozarimazetnthaaema . . . . ................ § ............ ,-"..
0,.“ ............... ................ ................ ............ -
0"p“..............é...”..-.........é................g...............% ................ :- ............ q
5 02 4.2 E i
. a...“ ................ ........................... _
I : : 0068: g
or....... ... ................§.............. .............E... ........§. .............. d
4.2.. ........................E ................................................. : ............... .1
4.‘l—.. . ................................................................ ..
4.6 ...............-................§ ................................................................. _.
4% l l. L l l
- 5 -1 —o.5 o as 1 15
I

Figure 4.9 Qualitative nonlinear function of the experimental system.

the identified values into the algebraic equations Aa = 4 formed from balancing the first

harmonic, we find that the force amplitude 23 to be 0.25 on average by using the first ten
extracted periodic orbits.

We also conduct a bifurcation diagram using the identified model by slowly increasing the
forcing amplitude as done experimentally, and sampling the steady-state response at the
same time interval. We carry out two bifurcation diagrams, one with the nonlinearity up to
the cubic term in the equation, the other one with the fourth term, as shown in Figure 4.10.
Both parameter sets came from the same 4th-degree model. Figure 4.10(a) shows a
nonperiodic response when the force amplitude is in the range of O.23~0.29 for the model
with cubic nonlinearity, and about 0.2~0.29 for the model with fourth power nonlinearity

in Figure 4.10(b). In the forcing range shown in the diagrams, the bifurcation diagrams are

72

 

 

 

 

 

 

 

 

‘ T Y '7 1 1 T
E 2 r ........................... |uuemzztt°.°.°.'.°. ....... '.'.’.'.'.‘.'.'.°.'.'.'.'.'.‘.’.'.'.°.‘.’.'.'.'.°.°.'.'.'.°.'.'.'.'.'.°.'.'77:.'.'.'.'.°.'.'.'L
1.00.00.003. !!:. '!.51
o .. .............................................................................................. ..
a I ......... :41 111i
§_2 .- ..................... ...? ..................................................................... ‘
4 1 l l l
0 0.1 0.2 0.3 0.4 0.5 0.8 0.7 0.8 0.0 1
W amplitude
‘ 1 17 j v— 1 T 7
g i5 ooooooooo 9 ooooooooooooooooooooooooooooo 6 ooooooooo O ooooooooo 9 ooooooooo 1
2.- ........ ', .................. .p....o....:. ......... . .......... :..........: .......... :.........‘. ........ ..4
. .1 :l; i": 3 s 3 s a s
‘ ------- o. : : : : : :
a o.,|:i ..i :II'::: ,. .......... .......... ... ........ _.
, .111 1111.”; 2 a s 2 z
I. . . . . . .
§_2. ......................... “1.! ......... -.........:. ......... .... ........ ..
_4 3m ; . t L : 1 i :
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
tom amplitude

Figure 4.10 Bifurcation diagram of the identified model, (a) with the cubic nonlinear
terminthe model, and(b) with thefourthpowerterminthe model

similar. As the forcing amplitude increases beyond four, the bifurcation diagrams are
different, since the fourth nonlinear term makes the system globally unstable.

The calculated forcing amplitude can generate chaotic motion, as happened
experimentally. Hence, a forcing amplitude 23 = 0.25 is used as a typical one in Eq. (4.11)
for numerical simulations.

Numerical integration of Eq. (4.11) is carried out using a Runge-Kutta method as usual.
The phase portrait of the model with cubic nonlinearity is shown in Figure 4.11, from
which the unstable periodic orbits are extracted as shown in Figure 4.12.

There is some resemblance between the periodic orbits extracted from the reconstructed
attractor in Figure 4.12 and the identified one in Figure 4.8. A large difference between

periodic orbits extracted from the experimental and the numerical models does not

73

 

Phuo space of tho identified model

r

 

 

.M

_
.

 

 

1.5

-0.5

-1

—2

0.5

—1 .5

 

5.5

 

1.5-

-1-

—1.5 -

Figure 4.11 Phase portrait of the identified model with cubic nonlinearity

 

pd-Z

 

 

 

 

 

—1

 

 

 

 

 

 

 

 

 

 

—‘I

Figure 4.12 Some periodic orbits extracted from the phase portrait of Figure 4.11

74
exclude the possibility of the existence of more similar orbits. Thus, we can conclude that

Eq. (4.11) qualitatively represents the experimental system. Given the mass and the scale

factor 7 discussed in section 4.6, the physical parameters can be obtained accordingly.

4.8 Estimation of the Natural Frequency and the Damping Ratio
Linearizing the identified model of Eq. (4.11) around the equilibrium points, we can
calculate the eigenvalues of the linearized model, and hence estimate the natural

frequencies and damping ratios for comparison with experimental measurements.

The Jacobian of Eq. (4.11) is

Df = 1 0 2 3 11 (4.12)
- 131- 2132x- 333x + 41341: -a

The equilibrium points are obtained from Figure 4.9 by locating the zero-crossing of the
nonlinear stifiness function. The equilibria for the curve with four nonlinear terms are
0.72 and -1.1. Then the eigenvalues of Eq. (4.12) are — 0.0175 2t 0.6348i and
— 0.0175 :1: 0.9426i in the time-normalized system. The real part represents the decaying
rate, and the imaginary part represents the undamped natural frequency (9.. . The damping
ratio can be calculated by dividing the real part by the imaginary part, yielding 2.76% and

1.86% for the right and the left well respectively. Converting to real time system by

multiplying by the driving frequency (7 .5 Hz in this case), the damped natural frequencies

rod = con 1 — c2 are 4.76 Hz and 7.07 Hz, respectively.
Omitting the fourth degree nonlinear term, the equilibrium points are 0.68 and -1.2 from

Figure 4.9, and the eigenvalues of the linearized model are —0.0175 i0.6541i and

75
-0.0175 :l:0.8626i . The damping ratio becomes 2.68% and 2.03%, and the real time

damped natural frequencies are 4.91 Hz and 6.47 Hz for the right and the left well
respectively
The transfer functions of the experimental system were measured by Bart Kimblel, as

shown in Figure 4.13.

 

m I I I I I I
a- ......... 6dB ......... 4

I C .

. . .

. .

l O

. O
. . .

. . .

. . .

. .

. .

b ooooooooooooooooooooooooooooooo
. . .
. . .
u u .
. u .
. . I
o o .
. . .
. . .
. o
.

.
.

I
l

    
      
       
  

-
I
I
I
I
I
I
O
I
I
I
I
I
I
I
I
Q
I
I
O
I
I
I
I
I
I
I
Q
I
I
I
I
I
I
I
I
I
I
Q
I
I
I
I
I
I
O
I
I
I
Q
I
I
I
I
I
I
I
I
I

l

O
I
i

 

 

 

5 ‘ -
5 .......;.l ........................................................... _
5 I
o ...-... ' ..........-........;..' ........................................................... q
_5 ......... g .......... I .........
1 : : : : :
-10... ......... . .......... .......... .1.. ........................................
5 : E I (0.; = 3.5112 . a -
-15!- ......... 3....."..IE.-....'...E...... ...E-ooco.o..-3.........II:.....-....E-on-o-ca..4:~.
: ' : | 3 : - : :
1 L L .1 L l l
202 4 6 8 10 14 16 18 20

Figure 4.13 Transfer function of the right well of the experimental system

Using the half-power'point method. and assuming the damping ratio, C, is small, the
damping ratio can be estimated by

(”Z-ml

20),.

C:

 

(4.13)

 

1. 'I'helrmsferfunctiomoftheexperirnerml systernwerekindly givenbme KimbleItPerm State University.

whet

76
where (a)1 and '33 are the frequencies at the half-power point of the transfer function, 0)..
is the undampec natural frequency, as indicated in Figure 4.13. Table 15 contains a
comparison of the natural frequencies and the damping ratios of the linearized model and
the experimental measurement. Since the stiffness function in the identified model, Eq.
(4.11), is in some sense a best fit to the experimental stiffness function, variation in the

slope of these functions leads to variations in linearized quantifies.

Table 15: Comparison of the natural frequency and the damping ratio

 

 

linearized linearized

nonlinearity nonlinearity
g 2.73% 2.68% 2.76%

left well “’4 8.5 Hz 6.47 Hz 7.07 Hz
g 2.52% 2.03% 1.86%

 

 

 

 

 

 

 

4.9 Discussion

Using the experimental data, we reconstructed the phase space of the experimental
system, from which the unstable periodic orbits were extracted for use in our parametric
identification algorithm. A qualitative model is obtained to represent the experimental
system as Eq. (4.11).

There are some discrepancies in the identification results when using different extracted

periodic orbits individually. We used several periodic orbits together because diflerent

For

M

77

periodic orbit VialiS different area of the reconstructed phase space. Using several of them
provides a better representation of the system behavior, and also gives some statistical
information of the identification results.

For this reason, chaotic dynamics is beneficial for parameter identification, particularly

when the form of the model is approximate. We will revisit this issue in next chapter.

CHAPTER 5

Errors in Parameter Estimates

5.1 Introduction

Errors are inevitable in any parametric identification method, arising from the incorrect
modelling, data acquisition and data manipulations. Modelling is not only a main source
of error, but also a critical factor for the success of an identification method. We have
discussed this issue when we presented our method in Chapter Two, and some criteria
have been used to validate the model in Chapter Three.

In this chapter, we will focus on the quality of the data, and its effect on the formulation of
the least-squares estimation Au = q . An obvious source of error is the noise, which is
inherent to the data acquisitions and manipulations. We will treat it as external to the
system response, and assume it is random and bounded. We have also used the unstable
periodic orbits exclusively in our identification scheme for a chaotic syStem. The unstable
periodic orbits are extracted from a chaotic attractor, and used as an approximation of the

real periodic orbit of the system. The deviation of the exuacted unstable periodic orbits

from the real one is another source oferror to be discussed in detail.

78

 

5.1

79
5.2 Errors induced by Noise

Suppose some noise Mr) is added to a periodic orbit p(t), such that

130) = 190) +110). (5.1)
where the noise is assumed to be uniformly distributed and uncorrelated to the periodic

«bit. Applying this noisy periodic-orbit data to a power-series model (assumed to be

valid),anonlineartermxk inthemodelisexpressedinaFourierseriesinaform

l

F{(.3)"} = F{p"+k..p*“... = 51/} +kF{np*‘ }..., (52)

where F denotes the Fourier series representation.
Suppose that the upper bounds of p(t) and n(r) are known, i.e., there exists real positive

numbers so and §,suchthat

Ip(t)l$5o, and,ln(t)lS§. (5.3)

for r e [0, T] . The Fourier coefficients of the real periodic orbit are
_ .2. f . . ..
ajk - 1' Op (1) cos (you)
2 (5.4)

k . .

bj,‘ = Til-p (t) srn Umt)dt
where the subscript it refers to the nonlinear term 1*.

Similarly, the Fourier coefficients of the noise-contaminated periodic orbit can be
calculated as

a}. = %1: 5‘0) cos (mod.
- 2 r (55)
bi. = 71:5 (t) sin (mod.

80
Thus, the difierences of the Fornier coefficients between 13‘ (t) and pk (I) will be

bounded by the following relationship:

Aajki alajk-erd = %11:(pk(t) —iik(t))cos (irot)dn1

k
2 T k! k-i i .
" 7110 Elm" m" (r) ”30“)“1 (5.6)
k o
It! It ‘
52:21 (t-i)!.‘!“’ ('25)

for every harmonic j. Normalizing the error bound of Eq. (5.6) by 50" yields

321.4...(g)....-.,(.g)2.... (.7,

The arguments will be the same for 111ij = |bj.— 5,4.

For It = 1 , the bound depends on the bound of noise only, which is usually assumed to be
as small as a few percent, specified by the noise-to-signal ratio, C/ p . As It increases, the
error bound will accumulate, making the perturbations in A larger accordingly. If the

degree of the nonlinear term in the mathematical model is excessively large, the accuracy

of the estimation results will deteriorate rapidly by the noise.

Having these bounds on 1Aa 1.11 and IAij , we can determine a bound on IIAAI and 1A4]

due to noisel. Then by a method given in swtion 5.4 below, a bound on the error of the
identified parameters, IAal , can be estimated.

 

1. AlluormeherearefiuclidemZ-norms.

Pl

81
5.3 Errors Induced by the Periodic Orbit Extraction

Anather source of errors comes from the extraction of the unstable periodic orbits. We

have specified a spatial criterion 6 in the state space, such that Ika-x‘ Se,tolook for

the recurrent points in a large chaotic data set. We presume that when the trajectory comes
close to a periodic orbit, it approximately follows the motion of that periodic orbit, so that
whentherecmrentpointsuelocatedthesegmentofdataistakenastheapproximate
periodic orbit. This approximation is related to a characteristic quantity of the associated
periodic orbit as discussed below.

Consider a neighborhood It (1) of a saddle fixed point 3 of a period-k orbit in the

Poincare section, as sketched in Figure 5.1.

     

Es : stable manifold

E“ : unstable manifold

 

 

Figure 5.1 Close look of the periodic orbit extraction on the Poincare section

82
Thedynamicsofthe systeminthisneighborhoodcanbeviewedintermsofalinearmap

T, such that

ka-it-an—x) , (5.8)
where T is the linearized map about the period-k orbit 2 in the period—k Poincare section.
which is invertible, since the periodic point is a saddle for a chaotic system, and x is the
saddle point, representing the true periodic orbit of a nonlinear system. If the orbit is near
the saddle to be considered as an approximate periodic orbit of period It, then by the

linearized map, the spatial distance between an orbit and its k iterations can be written as

Pk—xol = |(x,‘-1)- (xo-X)|
- [Toto—x) — (Io—1" (5.9)
a [(7—1) (xo-X)|Se

Taking the matrix norm, Eq. (5.9) is bounded by the singular values of the mauix (T - I) ,
such that

2“Mo-33' $I(T-’) (Io-1H Shh—3|. (5.10)
where k, and 3.2 are the maximum and the minimum singular values of the matrix
(T-l) . Since |(T-—l) (xo-x)| Se, the distance between the approximate periodic

orbit and the true periodic orbit, by the criterion Ix, - 10' S e , is bounded by

151 s lxo 4| 5 e/xz. (5.11)
Note that, since 35' is a saddle point, the singular value 2.2 will not be zero in any case.

Similarly, taking the map backward, the spatial distance becomes

83
I‘t‘xd = '09:") _T_!(xk-x)| = “fulfill-1”. (5.12)

and the distance between x,‘ and the saddle is bounded by

82: I‘r’l Se/uz. (5.13)
where rtz is the smallest singular value of the matrix (T4 —I), which is non-zero
becausei'isasaddlepoint. Thus theapproximateperiodicorbitbythecriterionof
l‘k‘xol Se willbe bounded by the largerof 81 and 82, i.e.

8 = max{81, 52} Smax{—e-, _e_} . (5.14)

MM

Since 12 and “2 are the minimum singular values of matrices relating to the linearized

map of the periodic orbit, the bound 8 in Eq. (5.14) bounds the error of the extracted
periodic orbit. The smaller the pro-specified spatial criterion 8 is, the smaller the bound 8
will be. However, the criterion 8 should be determined by the data set, since excessively
small 6 may result in no data points fitting in the criterion.

With the error bound 5 determined by the characteristic of the periodic orbit as above, we
canproceedtoboundtheerrorsinAandq,i.e. w-e .q'rfindaboundon [MI and IAqI
due to the periodic-orbit extraction. The argument is the same as that presented in previous
section, Eq. (5.2) to (5.6), except that the upper bound of noise C is replaced by the error
bound 5 .

In reality, we have to analyze the data to approximate the map 1'. The procedure follows
ideas of [5,18, 19,35,491.

 

 

 

Figure 5.2 A sketch of the construction of a linearized map

Let {xi} ,1 = l,...,P,beaset ofdatapointsinstate space, denofingaperiodic orbit

extracted from a chaotic attractorl. We intend to find a sequence of local linearized maps

Ti’j a 1, ...,M, by increments ofm, such that

xi”. =- zjj (5.15)
where M = P/m, an integer, and 10“,)“ = x", since the orbit is periodic. Then we

construct a compound map by multiplying each local linearized map T j in a reverse order,

such as

r = r,,rg,,._,...r1 (5.16)

This compound map T represents the torn! linearized map of the periodic orbit, as the one
used in Eq. (5.8) on the period-k Poincare section. Figure 5.2 is a sketch of the
consu'uction of this linearized map.

 

1. ItshouldbenoedthatweuaedtheumemtationxforpointsonthePoineIelectioneariierinthiuection.

85
Considera small ball ofradius rcenteredatthe orbital point xi, anda setofndata points

{xk'},i = l, ...,n, included inthe ball. The numberofdatapointsnis tobediscussed

later. The set of displacement vectors between 1*, and If is formed as

{y.-} = {xg- ill’h'lesr}' (5.17)

After the evolution of a time interval 1 a mm, the orbital point x]. will proceed to x,

' + Ill
and the neighboring points x“ to xk‘ + m . The displacement vectors y,- are thereby mapped

to

{1g} 3 {xk‘+m-xj+u|l.xk‘-x1|5r}. (5,13)

The evolution time interval mAt efi'ects the quality of the linearization, depending on the
dynamics of the system. If mAt is too small, the map will resemble an identity map. If
mm is too large, nearby points evolve beyond the regime of linearity. We have found no
proposed method for optimizing the choice of m. Usually we choose m to be one or two.

Let Y be a matrix containing the vectors y‘- and 2 containing the vectors 2‘. . Ifthe radius

r is small enough, evolution of Y to 2 can be approximated by a linear mapping

represented by a mauix Ti , such that

Z E T)! . (5.19)

By minimizing the squared error norm of Eq. (5.19) with respect to all components of the

matrix T}. , we obtain an expression for Ti as

86
I -l o
T.- = (W) (n). M»

where (.)' denotes the matrix transpose. The matrix T}. is an approximation of the
linearized map at xj. If there is no degeneracy, Eq. (5.20) can be solved for T1” In some

situations, Y may not be uniformly distributed in all direction of x]- (because of the data).

. 1
andmaynot spantheentirephase space, suchthat ()7). doesnotexistandthemauix
Ti may not be well defined. In such case, we might simply skip it and proceed to the next

point 1““, or use larger ball to find another set of neighboring points to eliminate the

degeneracy. In this way, we can construct a sequence of linearized maps for the periodic
orbits, and multiply each one around the periodic orbit in reverse order as that in Eq.(5.16)
to obtain the final map for use in Eq. (5.11) and (5.13).

To include enough neighboring points, n, around the orbital point for constructing the
linearized map, the radius of the neighborhood has to be larger than the spatial criterion
used to find the approximate periodic orbit. Lathop and Kostelich [35] used a radius of 68
to include 50 or more points. Eckmann and co-workers [5,18,19] increased the radius until
30 or more points were found, while Sano and Sawada [49] set the number of points to be
20, and confirmed that lower number still gave similar results, provided that number was

greater than the embedding dimension.

5.4 Sensitivity of the Parameter Estimates to Errors

Knowing the error bounds IAAI and IAqI induced by the noise and by the periodic orbit

extraction, we could procwd to estimate the sensitivity of the parameter estimation results

87
tothese errors.

Let AA and A4 betheperturbationsofAandqrespectivelyinthelinearsystem, withA

being a full-rank mauix. Define

a = qu
a’ = (A+AA)1(q+Aq) . (521)
r =- q-Aa

where A? is the pseudo-inverse of A, which yields the optimal least-squares solution of

Ar: = q,asdoes the (A +AA)T for a“ inEq. (5.21). Assuming

AAI IAqlz ___l__
"'"m [A ’25] ]< (—T)cand (512)
and
Id;
° 9 <1, (5.23)
sm( )ilzl;

implicitly defining 9, 0 $0<rt/2, then theerrorin the parameters a is bounded by [4]

(1* — a
I '2 Sn [M + tane [cond (A) ] 2] + 0(112). (524)
2 cosO

For a given model and given periodic orbits, A and q are fixed. and ‘n is determined by
IAAI and IAqI induced by the error sources discussed in previous sections. If the model
is chosen properly, the residual, r, is usually very small, resulting in small 9 . Thus the
error bound is Eq. (5.24) depends linearly on coud(A). If the model is not properly chosen,

the error bound is Eq. (5.24) will depend on the square of cond(A), making the error bound

sufficiently larger. From numerical results presented in Chapter Three, we found that the

88

choice of a model is the most influential factor on estimating the parameters. Also, an
improperly chosen model results in a large condition number for A, resulting in a violation
of the assumption in Eq. (5.22). In fact, this occurred when we tried to compute the bound
of the error in section 3.1.5. In that example, we directly computed IAAI and Mg] by
comparing the quantities in the noise-free and noisy cases. Thus, the estimate (5.23) may

not always be practical, but at least indicates some trends.

5.5 Using Several Periodic Orbits

In this section, we discuss how an application of several periodic orbits together in the

identification scheme, as opposed to using a single periodic orbit, can reduce the

sensitivity of the parameter estimates to errors.

Suppose a mathematical model has been chosen, and difi‘erent matrices A‘. are formed

using p different periodic orbits, such that i = 1,...,p. Suppose different sets of
parameters (“1,- are then estimated from
Aiai=qi’i = 1,...,p. (5.25)

Then, combining several A ‘- into a single matrix A, such that

=[A1,A2,,A1], (5.26)
and combining the corresponding q‘. into a single vector q, we obtain a set of parameters

(‘1 from the combined equation, At“: = q , through

LATAJG = ATq, (527)

which is equivalent to

89
P P
T T
[ZAJ'AI]& = 2AM
j-l j-r

Summing Eq. (5.25) on i, and subtracting from Eq. (5.28), yields

iAfAj(a—&j) = £11,7er (at—a) +(a-&j)) = o.

i-1 i-1

where a is the true parameter vector. Defining

we obtain an expression of the error using different sets of periodic orbits, as

P

A :3 le-Aj
j-l

(5.28)

(5.29)

(5.30)

(5.31)

9
Note that, z R]. = 1. Thus, Eq. (5.31) implies that the error in parameter estimates from

1'1

several periodic orbits is a weighted average of parameter errors from individual periodic

orbits, through the weights R j. If a particular periodic orbit were known to yield the best

estimation, we would use it for identification. However, we have no idea which one will

be the best in general. Hence, using a combination of several periodic orbits has

advantages. It can give a reasonable estimate, reduce the sensitivity to errors, and improve

the statistical properties of the identified results, which has been shown in our previous

applications.

CHAPTER 6

Conclusions and Future Work

6.1 Conclusions

We have presented a method for identifying parameters of nonlinear systems that exhibit
chmdcbehavimthisismexmsimofmexisfingmerhodfunodinursymswim
mbkpaiodiemponnWeexploitthechwdcamcmrofthenmfineusystemby
exuacfingmeunstablepefiodicabiufiomthechmficamacmrmmpreaentthesyswm
behavimfiachtammthemamemaficalmoddisexpmssedmafimmFomiersaiesusing
meexuacwdpaiodicorbinaudtheharmonic-bdancememwuappfiedmfamasuof
linear algebraic equations in system parameters for least-squares estimation. We have
demonstrated that the present approach is applicable to externally excited. parametrically
excited. and autonomous chaotic systems. This may not be feasible using Other methods.

Alfioughchmshasbeenmgadeduundesimbkmisembediscardedfiomthephysical
systemsithasrichinformationcontentascomparedtoaperiodicuajectory.Thisrichness
hasbeenexplouedindimensionafitysmdies,nmhneupredicdomandcmuokm
cenualmemeismepresenceofunsmbkpuiodicmbimwhichcmbeexmwdandused
mchamctuizemechaoficamacmr.1haefmemeymusefidinpanmeuicidendficadon
faachmdcsysmthammetbeavaflabihtyofmanypaiodicabitswimmfiuem

90

91
puiodsprovidesa‘pasismntexdndon’fuasyswmwithmanypuammsism
advantageofusingchaosinparametricidentification.
Modemngismimpammisminmepanmuicidmfificafimmmmmust
beablemeamhmeeuendflchuactuisdcsofmesysmundainveefigafimwhume
famofthenmfineafityandthetypeofexeimfionFamoothsyswmsthroughommis
smdy,wechooseapowasaiermnpresentmesymnmfinafity,basedmmefacttha
aMfimcfimcanbeexprusedinapowasaieemndtheequifibnumpoinuSome
quecdmsrdatedmthisreprecenmfiomsuchuhowmcboosetheopdmaluuncadmof
mepowasaieaandwhethaunmthedanuewiminthendiusofemvagmceofthe
“remainmbesmdiedfmmezWehawdiscussedmispmblemusingnnumuical
exampleinChapterTwo.Wefoundthat.whenthedataarewithinthendiusof
convagenceofthepowasaieamemodelcanbeaccmulyideudfiedWhenthedanm
omofmeradiuofwnvagemednmmyofidenfifiadmreculudewrimhsuch
usethemdnonfineufirnefimofthesystemhasmbeknowninadermobninmte
results. Knowing something about the system under investigation is fundamental for
pmmeuicidenfificafionOtberwisewehavembecmwntwimpufialdesuipfionofthe
identified system. The moral is that good models lead to good quantitative results.
Wehaveappfiedmismethodnunnfieanymseveralchaoficsystemasuchuafmced
Dumngoscinamr,asmothoulomb&icdonoscmamr,apanmeuicanyexciwdbcam.
andalaenzoscillatoe‘l'heacctuacyistypicallywithin 1%erroroftheidentified
puametusfmnoise-fieedehenmemiselwelinaeasamewcmacydemetes,
although in arobust way, especially when the mathematical model is not properly chosen.

Wehavealsoappfiedthismethodtoafacedmechanicaloscinamrwithatwo-wefl

92
sfifi‘nesspotenfiaLChmficdanwaewquiredfiomasuaingaugebymmbg-todigiml
convahedamwuepmpwfioualmthedisphcementofmeosdlhtmlheamphmde
oftbepuiodicfmcingfuncdmwuunknownWemdeleddteexpaimennlsymina
powasaiealheunstablepaiodicubiuwaeexuacwdfiomareconsuuctedm
usingmemethodofdehys.Cmsimcyofmeidenfifiedpmmrwuachievedusing
difierentcombinafiomofmeexmcwdpaiodieabimladingmaquafimfiwduuipfion
ofmesynemnmfineafity,andhencemaqualimfiwmoddoftheexpuimenmlsymA
bifineafimdiagramwucmsuucwdmingtheidenufiedmodeaneafingphenomenaof
period doublingandchaos. 'l'hisqualitativemodelcanbeusefulforfru'therinvectigation
oftheexperimentalsystem.
ModekwaevaifiedbycompuingtheLyapumvexpooenmbifmcadondiamand
mesuucuneofthemsmbkpaiodicabiminmemigindandtheidenfifiedsymWe
almsoughtconsistencyinmeidendfiedpanmembyusingdifluemsenofpaiodic
orbits.
Nomainsomcesofamrsinthepanmetaesfimateamiseandtheexuacfionofthe
unstable periodic orbits, have been examined closely from statistical and geometrical
pointsofview.Weconsideredthenoisetobeuncorrelatedtotbesystemrecponseand
uniformly distributed. We constructed a linearized map for the periodic white, and derived
aboundbythe singularvaluesofthecorrespondinglinearizedmap.Byexpressingthe
noisyperiodicorbitsinFourierseries,theerrorintheFomiercoeficientsofanonlinear
taminthemamemaficdmodelcanbeboundedasafumfimofthemise-m-signalnfio.
Emencrmealgebraicequafiominammmeufashiomhmreasingdwuncanintyof

meideufificafionrendmfihemboundofthepmmetaesdmtesispmpadmalmthe

93
squareofthecondifionalnumberofthemauionftheidenfificafionequations, Aa - q,

ifmemsidualsamhrge.0mawise,theurorboundishneuincond(A).Wefmmdthm
udngsevaalexuacwdpaiodicabinmgemawinreducemesensifivityofmepmmem

ertimatestotheeeerrors.

6.2 FutureWork
Wehaveusedmemsmblepuiodicubinexdusivdywrepruentthesyscmbehaviorm
mnpanmeuicidenfificafionschemeJhesynemisinafumofmadimydifl'aendal
equation,andrecastasasetoflinearalgebraicequafionsthroughthebalancingofthe
Fomiercoeficientsofeachtermintheequadon'l‘hisisstrategicallyconvenientfa'our
appficafimabmnmmathemadcanyfigorouthntha’smdiermyleadmabeuer
understandingofhowitrelawatotheordinarydifierentialequation.
WeexpmssedthepuiodicmbitinaFomiasaiegudngmefundmnnlfiequencyin
cdmhfingtheirFomiacoefficieutalheaniersaiecofapaiodNabitcmMof
frequency components m/N, m = 1,...,N, thus making available manyharmonicsto
balmceJtmaybewonhexamimgwhethasubhumonicsmmpahumomcsyieldbeuer
results.

Errors in the identification process can be investigated more thoroughly by examining
howthenyestmemselvuinthefmmulafionofthelean-aqumesfimafiomswhu
theexnacfionoftheunsmbleabiufiommechaodcamacfingsetandthefimmafionof
the Foruier coefficients of these unstable periodic orbits.
Wechosetousepolynomialsasthebasisfuncfionstomodelthesystemnonlinearity.
Other form of basis functions such as wavelets, Pade functions, sigmoid functions. and
mdialfuncnonamightbeworthinvestigation.

94
Expaimmmlworkcanbeexplaedfmmer,suchuwithmmnomoussymmsandhigha

degree-of-freedom systems. Quantitative studies, involving specific material or system
pmmewrxmhuinidenfifyingmeduficmoduluxnfightbevduabkExpahneoml
audiumdisuibuwdsyswmswouldcomplemntthenumaicalwakonuudaandco-
workersl59, 60, 61],andmightraiseotherinteresting'issues.
ThissmdyfocusedonsmoothnmfinarfumdonaNm-moothsystemasmhuimpact
andfiicfiommightcallforadjustmentsinthismethod.
Ammpfingmkmemodmidwfifypmmamwhichmexpressednonlineuiyinme
difi‘erential equation of motion, such as in arguments of functions (e.g. tanh(a.x))
requimssomedevebpmenLAnexcusimofthepresentmahodmmoadbandey

balandngtheredandimaginaryparuofdreFomiermsfumofeachmmay
accommodatethecaseofnoisyinput.

BIBLIOGRAPHY

[1]

[2]

[3]

[4]

[5]

[5]

[7]

[3]

[9]

95

BIBLIOGRAPHY

Abarbanel, H. D. 1., Brown, R., and Kadtke, J. B., “Prediction in Chaotic Nonlinear
Systems: Methods for time series with Broadband Fourier Spectra”, Physical Review
A, 41, 1990, 1782-1807

Abarbanel, H. D. 1., Brown, R., and Tsimring, L. S. “The Analysis of Observed
Chaotic Data in Physical Systems”, Reviews of Modern Physics 65, 1993, 1331-
1392

Aguirre, L. A. and Billings, S. A., “Validating Identified Nonlinear Models with
Chaotic Dynamics”, International J. Bifurcation and Chaos, Vol. 4, 1994, 109-125

Atkinson, K. E., An Introduction to Numerical Analysis, 2nd ed., John Wiley and
Sons, 1989

Auerbach, D., Cvitanovic, P., Eckrnann, J. P., Gunaratne, 6., and Procaccia, 1.,
“Exploring chaotic motion through periodic orbits”, Physical Review Letter, 58,
1987. 2387-2389

Beck, J. and Arnold, K., Parameter Identification in Engineering and Science, John
Wiley & Sons, 1977

Breeden, J. L., and Hubler, A., “Reconstructing Equations of Motion from
Experimental Data with Hidden Variables”, Technical Report CCSR-90-7,
University of Illinois at Urbana-Champaign, 1990

Bridge, J. and Rand, R., “Chaos and Symbol Sequence in System with a Periodically
Disappearing Figure-Eight Separatrix,” Biftu'cation Phenomena and Chaos in
Thermal Convection, ASME HTD-Vol.214/AMD-Vol.138, 1992, 47-55

Broomhead, D. S., and King, G. P., “Extracting Qualitative Dynamics fiorn
Experimental Data”, Physica D 20, 1986, 217-236

[10] Coppola, V. T., and Rand, R. H., “Chaos in a System with a Parametrically

Disappearing Separatrix,” Nonlinear Dynamics 1, 1990, 401-420

[11] Crutchfield, J. P. and McBamara, B. 8., “Equations of Motion from a Data Series,”

Complex System 1, 1991, 417-452

[12] Cusumano, J. P. and B. W. Kimble, “Experimental Observation of Basins of

Attraction and Homoclinic Bifurcation in a Magneto-Mechanical Oscillator,” in
Nonlinearity and Chaos in Engineering Dynamics, edited by J. M. T. Thompson and
S. R. Bishop, John Wiley & Sons Ltd, 1994

96

[13] Cusumano, J. P., and Sharkady, M. T., “An Experimental Study of Bifurcation.
Chaos, and Dimensionality in a System Forced Through a Bifurcation Parameter” To
appear

[14] Cvitanovic, P, Gunarame, G. H., and Procaccia, I., “Topological and metric
properties of Henon-type strange attractor,” Physical Review A, 38, 1988, 1503-
1507

[15] Ding, M., Grebogi, C., Ott, E., Sauer, T., and Yorke, J. A., “Plateau Onset for
Correlation Dimension: When Does It Occur?”, Physical Review Letters, 70, 1993,
3872-3875

[16] Ditto, W. L., Rauseo, S. N., and Spano, M. L., “Experimental Control of Chaos”,
Physical Review Letter, 65, 1990, 3211-3214

[17] Distefano, N., and Rath, A., “System identification in nonlinear structural seismic
dynamics”, Comput. Mech. Appl. Mech. Engng 5, 1975, 353-372

[18] Eckmann, J.,-P., “Ergodic theory of chaos and strange attractors”, Reviews of
Modern Physics. 57. 1985, 617-656

[19] Eckmann, J.,-P., Kamphorst, S. O., Ruelle, D. and Gliberto, S., “Lyapunov
Exponents from Time Series”, Physical Reviews A, 34, 1986, 4971-4979

[20] Eisenhammer, T., Hubler, A., Packard, N., and Kelso, J. A. 8., “Modeling
Experimental Time Series with Ordinary Differential Equations”, Technical Report
CCSR-89-7, University of Illinois at Urbana-Champaign. 1989

[21] Farmer, J. D. and Sidorowich, J. J., “Predicting ChaOtic Time Series”, Physical
Review Letters, 59, 1987, 845-848

[22] Feeny, B. F. and Liang, J -W, “Phase-space reconstruction of stick-slip systems”,
submitted to the 15th Biennial ASME conference on Noise and Vibration, 1995

[23] Feeny, B. F., N ayfeh, A. H., and Mook, D. T., “Parametric identification via the
method of multiple scales: modeling a ship at sea”, to appear

[24] Feeny, B. F., and Moon. F. C., “Chaos in a Forced Dry-Friction Oscillator:
Experiments and Numerical Modelling”, J. Sound and Vibration, 170, 1994, 303-323

[25] Fraser, A. M., and Swinney, H. L., “Independent Coordinates for Strange Attractors
from Mutual Information”, Physical Review A, 33, 1986, 1134-1140

[26] Gershenfeld, N., “An Experimentalist’s Introduction to the Observation of

Dynamical Systems”, in Directions in Chaos, Vol. 2, ediwd by Hao B.-L., World
Scientific, N. J ., 1988

[27] Gottlieb, O., Feldman, M., and Yim, S. C. 8., “Parametric Identification of Nonlinear

97
Ocean Mooring Systems Using the Hilbert Transform”, preprint

[28] Grassberger, P. and Procaccia, 1., “Characterization of Strange Attractors”, Physical
Review Letter, 50, 1983, 346-349

[29] Grebogi, C., Hammel, S. M., Yorke, J. A., and Sauer, T., “Shadowing of Physical
Trajectories in Chaotic Dynamics: Containment and Refinement,” Physical Review
Letters, 65, 1990, 1527-1530

[30] Guckenheimer, J ., Holmes, P., Nonlinear Oscillations. Dynamical Systems, and
Bifurcations of Vector Fields, Springer-Verlag New York, 1983

[31] Hammel, S. and Heagy, J., “Chaotic System Identification Using Linked Periodic

Orbits”, SIAM Conference on Applications of Dynamical Systems, Snowbird. Utah,
1992.

[32] Hanagud, S. V., Meyyappa, M., and Craig, J. 1., “Method of Multiple Scales and

Identification of Nonlinear Structural Dynamic Systems” AIAA Journal, Vol. 23,
May 1985, 353-372

[33] Ibanez, P., “Identification of Dynamic Parameters of Linear and Nonlinear Structural
Models from Experirrrental Data,” Nuclear Engineering and Design, 25, 1973, 30-41

[34] Kennel, M. B., Brown, R, and Abarbanel, H. D. 1., “Determining Embedding
Dimension for Phase Space Reconsu'uction Using a Geometrical Construction”,
Physical Review A, 45, 1992, 3403-3411

[35] Kesaraju, R. V., and Noah, S. T., “Characterization and Detection of Parameter

Variations of Nonlinear Mechanical Systems”, Nonlinear Dynamics 6, 1994, 433-
457

[36] Lathrop, D. P. and Kostelich, E. J., “Characterization of an experimental strange
attractor by periodic orbits”, Physical Review A 40, 1989, 4028-4031

[37] Ljung, L, System Identification-Jimmy for the User, Prentice-Hall, 1987
[38] Malasoma, J-M., Lamarque, C.-H., and Jezequel, L., “Chaotic Behavior of a

Parametrically Excited Nonlinear Mechanical System,” Nonlinear Dynamics 5,
1994, 153-160

[39] Mickens, R. E., “Bounds on the Fourier Coefficients for the Periodic Solutions of
Nonlinear Oscillator Equations”, J. Sound and Vibration, 124, 1988, 199-203

[40] Mohammad, K. S., Worden, K., and Tomlinson, G. R., “Direct parameter Estimation
for Linear and Non-linear Structures”, J. Sound and Vibration, 152, 1992, 471-499

[41] Mook, D. J. and Junkins, J. L., “Minimum Model Error Estimation for Pooriy
Modeled Dynamic Systems” Journal of Guidance, Vol. 11, No. 3, 1988, 256-261

98

[42] Mook, D. J ., “Estimation and Identification of Nonlinear Dynamic Systems” AIAA
Journal, Vol. 27, July 1989, 968-974

[43] Moon, F. C., Chaotic and Fractal Dynamics-An Introduction for Applied Scientists
and Engineers, John Wiley and Sons, New York, 1992.

[44] Moon, F. C. and Holmes, P., “A Magnetoelastic Strange Attractor,” J. Sound and
Vibration, 65, 1979, 275-296

[45] Nayfeh, A. H., “Parametric Identification of Nonlinear Dynamic Systems”.
Computer and Structures Vol. 20, No. 1-3, 1985, 487-493

[46] Nayfeh, A. H. and Mook, D. T., Nonlinear Oscillations, 1979

[47] Ott, E., Grebogi, C., and Yorke, J. A., “Controlling Chaos”, Physical Review Letters.
64. 1990. 1196-1199

[48] Packard, N. H., Crutchfield, J. P., Farmer, J. D. and Shaw, R. 8., “Geometry from a
Time Series”, Physical Review Letters, 45, 1980, 712

[49] Sano, M., and Sawada, Y., “Measurement of the Lyapunov Specu'um from a Chaotic
Tome Series”, Physical Review Letters, 55, 1985, 1082-1085

[50] Shinbrot, T., Ott, E., Grebogi, C., and Yorke, J. A., “Using Chaos to Direct
Trajectories to Targets”, Physical Review Letters, 65, 1990, 3215-3218

[51] Strang, 6., Linear Algebra and Its Applications, 1976

[52] Stry, G. I. and Mook, D. J., “An Experimental Study of Nonlinear Dynamic System
Identification” Nonlinear Dynamics, 3, 1992, 1-11

[53] Takens, F., in Dynamic System and Turbulence, Warwick 1980, edited by D. Rand
and L.-S. Young, Lecture Notes in Mathematics No. 898 (Springer, Berlin, 1981)

[54] Tufillaro, N. B., Abett, T., and Reilly, J., An Experimental Approach to Nonlinear
Dynamics and Chaos, Addison-Wesley Publishing Co, 1992

[55] Ueda, Y., “Randomly Transitional Phenomenon in the System Governed by
Duffing’s Equation,” J. Statistical Physics 20, 1979, 181-196

[56] Ueda, Y., The Road to Chaos, Aerial Press, 1994

[57] Wigdorowitz, B. and Petrick, M. H., “Modelling Concepts Arising fi'om an
Investigation into a Chaoric system,” Math. Comput. Modelling 15, 1991, 1-16

[58] Wolf, A., Swift, J. B., Swinney, H. L., and Vasano, J. A., “Determining Lyapunov
Exponents from a Time Series,” Physica 16D, 1985, 285-317

99

[59] Yasuda, K., Kawamura, S., and Watanabe, K., “Identification of Nonlinear Multi-
Degree-of-Freedom Systems (Presentation of an Identification Technique)”, J SME
International Journal, Series 111, Vol. 31, 1988, 8-14

[60] Yasuda, K., Kawamura, S., and Watanabe, K., “Identification of Nonlinear Multi-
Degree-of-Freedom Systems (Identification Under Noisy Measurements ”, J SME
International Journal, Series III, Vol. 31, 1988, 302-309

[61] Yasuda, K., and Kamiya, K.. “Identification of a Beam (Proposition of an

Identification Technique)”, J SME International Journal, Series 111, Vol. 33, 1990,
535-540

 

: lltl‘llu‘i: .

it). Ill. .

HICHIGQN S

TRTE UNIV. LIBRRRIES
111|111111111111111111111HI
9301 41 72708

1|“ l

 

312