Close
About
FAQ
Home
Collections
Login
USC Login
Register
0
Selected
Invert selection
Deselect all
Deselect all
Click here to refresh results
Click here to refresh results
USC
/
Digital Library
/
University of Southern California Dissertations and Theses
/
Multiple degree of freedom inverted pendulum dynamics: modeling, computation, and experimentation
(USC Thesis Other)
Multiple degree of freedom inverted pendulum dynamics: modeling, computation, and experimentation
PDF
Download
Share
Open document
Flip pages
Contact Us
Contact Us
Copy asset link
Request this asset
Transcript (if available)
Content
MULTIPLE DEGREE OF FREEDOM INVERTED PENDULUM DYNAMICS: MODELING, COMPUTATION, AND EXPERIMENTATION by Cheng-Yuan Jerry Chen A Dissertation Presented to the FACULTY OF THE GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Ful llment of the Requirements for the Degree DOCTOR OF PHILOSOPHY (AEROSPACE AND MECHANICAL ENGINEERING) May 2009 Copyright 2009 Cheng-Yuan Jerry Chen ii Table of Contents List of Figures iv Abstract x Chapter 1 Inverted Pendulum Dynamics 1 1.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Motivation and Related Early Work on Inverted Pendulum Dynamics . . . 4 Chapter 2 Inverted Pendulum Under Periodic Vertical Forcing 7 2.1 Single Degree of Freedom Inverted Pendulum Model . . . . . . . . . . . . . 7 Equations of Motion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 Linearized sDOF Model and Mathieus Equation . . . . . . . . . . . . . . . 13 Numerical Integration of the Nonlinear sDOF Inverted Pendulum Model . . 28 Nonlinear Analysis of sDOF Model Using Perturbation Theory . . . . . . . 41 Chapter 3 Experiments on the sDOF Inverted Pendulum 57 3.1 Material and Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57 Driving Inverted Pendulum by Loud Speaker . . . . . . . . . . . . . . . . . 57 Inverted Pendulum Motion Using High Speed Camera . . . . . . . . . . . . 60 System Calibration and Image Processing by NI Vision Builder . . . . . . . 63 Data Analysis and Data Filtering . . . . . . . . . . . . . . . . . . . . . . . . 64 3.2 Single Degree of Freedom System Experimental Results . . . . . . . . . . . 67 Chapter 4 Multiple Degree of Freedom Inverted Pendulum Model 74 4.1 Formulate nDOF System Equation of Motion From Lagrangian Approach . 74 4.2 Validate the Single Degree of Freedom Inverted Pendulum Model . . . . . . 81 4.3 Two Degree of Freedom Inverted Pendulum Model . . . . . . . . . . . . . . 82 4.4 Coe¢ cients of 2DOF Inverted Pendulum Equation of Motion . . . . . . . . 86 4.5 Numerical Integration on Nonlinear 2DOF Inverted Pendulum Model. . . . 89 4.6 Small Angle Assumption to Obtain Linear 2DOF Inverted Pendulum Model 92 iii Simulation Results on Linear 2DOF Inverted Pendulum Model under Mode 1 Vibration . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 Simulation Results on Linear 2DOF Inverted Pendulum Model under Mode 2 Vibration . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101 Numerical simulation result comparison . . . . . . . . . . . . . . . . . . . . 110 Chapter 5 Experiments on the 2DOF Inverted Pendulum 115 5.1 Material and Methods for 2DOF Experimentations . . . . . . . . . . . . . . 115 5.2 2DOF System Experimental Results . . . . . . . . . . . . . . . . . . . . . . 120 Experimental Results on 2DOF Mode 1 Vibration . . . . . . . . . . . . . . 120 Experimental Results on 2DOF Mode 2 Vibration . . . . . . . . . . . . . . 124 Experimental Results on 2DOF Upright Vertical Mode . . . . . . . . . . . . 128 Chapter 6 Summary and Future Work 133 Bibliography 135 Appendix 140 iv List of Figures 2.1 Single degree of freedom inverted pendulum under vertical periodic forcing at its base P: . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 2.2 Transition curves for the linear Mathieu equation on the plane separate regionsofboundedsolution(shadedarea)fromregionsofunboundedsolution or unstable system(white area). . . . . . . . . . . . . . . . . . . . . . . . . . 23 2.3 Stutts diagram for stability regions separated by the transition curves of n = 0 and n = 1 within the range of0:5< < 0:5. . . . . . . . . . . . . . 25 2.4 Stabilityregimediagramofsingledegreeoffreedominvertedpendulumunder vertical periodic forcing showing in the ! plane. . . . . . . . . . . . . . 27 2.5 Numerical integration results of the single degree of freedom inverted pendu- lum under normalized forcing amplitude = 0:1708 and normalized forcing frequency ! = 8:7313 with initial condition 0 = 0:1 and _ 0 = 0. (a)System time response. (b)Angular velocity time trace. (c)System phase portrait. (d)3D view of phase portrait with respect to periodic forcing. . . . . . . . . 32 2.6 Numerical integration results of single degree of freedom inverted pendulum under = 0:1710 and ! = 14:2784 with initial condition 0 = 0:1 and _ 0 = 0. (a)System time response. (b)Angular velocity time trace. (c)System phase portrait. (d)3D view of phase portrait with respect to periodic forcing. . . . 34 2.7 Numerical integration results of single degree of freedom inverted pendulum under normalized forcing amplitude = 0:2193 and normalized forcing fre- quency ! = 14:2784 with initial condition 0 = 0:1 and _ 0 = 0. (a)System time response. (b)Angular velocity time trace. (c)System phase portrait. (d)3D view of phase portrait with respect to periodic forcing. . . . . . . . . 35 v 2.8 Numerical integration results of single degree of freedom inverted pendulum under normalized forcing amplitude = 0:400 and normalized forcing fre- quency ! = 12:8506 with initial condition 0 = 0:56 and _ 0 = 0. (a)System time response. (b)Angular velocity time trace. (c)System phase portrait. (d)3D view of phase portrait with respect to periodic forcing. . . . . . . . . 37 2.9 Numerical integration results of single degree of freedom inverted pendulum under normalized forcing amplitude = 0:200 and normalized forcing fre- quency ! = 14:2784 with initial condition 0 = 0:1 and _ 0 = 0. (a)System time response. (b)Angular velocity time trace. (c)System phase portrait. (d)3D view of phase portrait with respect to periodic forcing. . . . . . . . . 38 2.10 Numerical simulation results of periodic solutions from varying normalized forcing amplitude with xed normalized forcing frequency ! = 14:2784. Graphs from top to bottom show values increase from = 0:1019 to = 0:1240. From left to right showing four geometric representations of system response: time trace, phase portrait, 3D phase portrait, and PSD, respectively. 42 2.11 (Cont. of 2.10) Numerical simulation results of periodic solutions from vary- ing normalized forcing amplitude with xed normalized forcing frequency ! = 14:2784. Graphs from top to bottom show values increase from = 0:1293 to = 0:1531. From left to right showing four geometric repre- sentations of system response: time trace, phase portrait, 3D phase portrait, and PSD, respectively. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 2.12 (Cont. of 2.11) Numerical simulation results of periodic solutions from vary- ing normalized forcing amplitude with xed normalized forcing frequency ! = 14:2784. Graphs from top to bottom show values increase from = 0:1611 to = 0:3346. From left to right showing four geometric repre- sentations of system response: time trace, phase portrait, 3D phase portrait, and PSD, respectively. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 2.13 Numerical simulation based on perturbation analysis solutions. All simula- tion is under xed normalized driving amplitude = 0:1233. . . . . . . . . . 54 2.14 Quantitativecomparisonoftheperturbationanalysissimulationresulttothe experimental result. (a) ltered experimental time trace and phase portrait (b)simulation time trace and phase portrait. . . . . . . . . . . . . . . . . . . 56 3.1 Single degree of freedom inverted pendulum experimental setup. . . . . . . 59 3.2 Preciseexperimentalmeasurementsfromimageacquisitionandimageprocess- ing using XS-3 high-speed camera, NI Vision Builder, MS Excel and Matlab. 60 vi 3.3 Three dimensional view of the SDOF experimental setup. . . . . . . . . . . 62 3.4 NIVisionBulidersoftwareappliedtoasingledegreeoffreedomexperimental setup. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 3.5 Zoom-in view of experimental data of inverted pendulum base motion. X denotes vertical periodic forcing displacement which has coordinate u(t) = "cos!t; Y denotes horiontal motion. . . . . . . . . . . . . . . . . . . . . . . 66 3.6 Experiment on single degree of freedom pendulum under 10:0 0:2mm driving amplitude and 25Hz driving frequency. System has normalized am- plitude " = 0:140:01 and normalized frequency ! = 13:30:2. . . . . . 69 3.7 Experimentonsingledegreeoffreedompendulumunder8:00:3mmdriving amplitude and 25Hz driving frequency. System has normalized amplitude " = 0:110:01 and normalized frequency ! = 13:30:2. . . . . . . . . . . 71 3.8 Experimentonsingledegreeoffreedompendulumunder6:50:3mmdriving amplitude and 30Hz driving frequency. System has normalized amplitude " = 0:090:01 and normalized frequency ! = 15:90:2. . . . . . . . . . . 73 4.1 Geometric relationship of gereral case nDOF inverted pendulum. under ver- tical oscillatory forcing at its base. . . . . . . . . . . . . . . . . . . . . . . . 76 4.2 Simpli ed two degree of freedom inverted pendulum model with two point masses mounted on two rigid rods. Two modes of vibration shown on graph: Mode 1(left) and Mode 2(right). . . . . . . . . . . . . . . . . . . . . . . . . 88 4.3 2DOFinvertedpendulumundermode1vibrationwithverticalperiodicforc- ing amplitude " = 0:00593 and forcing frequency ! = 250rad=s. Initial conditions are 1 0 = 0:05; _ 1 0 = 0 and 2 0 = 0:1; _ 2 0 = 0 . . . . . . . . . . . . 91 4.4 When applying di¤erent initial conditions to the same 2DOF inverted pen- dulum model we can get mode 2 vibration with the same vertical periodic forcing amplitude " = 0:00593 and forcing frequency ! = 250rad=s. Initial conditions are 1 0 = 0:1; _ 1 0 = 0 and 2 0 =0:101; _ 2 0 = 0. . . . . . . . . . 93 4.5 2DOFinvertedpendulumundermode1vibration. Theperiodicinputforcing amplitude was " = 0:2138 and the forcing frequency was ! = 200. . . . . . . 99 4.6 2DOFinvertedpendulumundermode1vibration. Theperiodicinputforcing amplitude was " = 0:03025 and the forcing frequency was ! = 200. . . . . . 100 vii 4.7 2DOFinvertedpendulumundermode1vibration. Theperiodicinputforcing amplitude was " = 0:03023 and the forcing frequency was ! = 200. The inverted pendulum shows quasiperiodic response. . . . . . . . . . . . . . . . 102 4.8 2DOFinvertedpendulumundermode2vibration. Theperiodicinputforcing amplitude was " = 0:2138 and the forcing frequency was ! = 200. . . . . . . 104 4.9 2DOFinvertedpendulumundermode2vibration. Theperiodicinputforcing amplitude was " = 0:1514 and the forcing frequency was ! = 200. . . . . . . 105 4.10 2DOFinvertedpendulumundermode2vibration. Theperiodicinputforcing amplitude was " = 0:11637 and the forcing frequency was ! = 200. . . . . . 106 4.11 2DOFinvertedpendulumundermode2vibration. Theperiodicinputforcing amplitude was " = 0:05037 and the forcing frequency was ! = 200. . . . . . 108 4.12 2DOFinvertedpendulumundermode2vibration. Theperiodicinputforcing amplitudewas" = 0:05andtheforcingfrequencywas! = 200. Theinverted pendulum shows quasiperiodic response. . . . . . . . . . . . . . . . . . . . . 109 4.13 2DOFinvertedpendulumundermode2vibration. Theperiodicinputforcing amplitudewas" = 0:26andtheforcingfrequencywas! = 200. Theinverted pendulum shows near periodic response. . . . . . . . . . . . . . . . . . . . . 111 4.14 Phaseportraitcomparisonofperiodicsolutionsonnumericalsimulationusing linear 2DOF inverted pendulum model. From upper left to lower right, the forcing amplitude has ascending values which result in descending frequency ratio. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112 4.15 Powerspectraldensitycomparisonofperiodicsolutionsonnumericalsimula- tion using linear 2DOF inverted pendulum model. From upper left to lower right, the forcing amplitude has ascending values which result in descending frequency ratio. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114 5.1 Two degree of freedom inverted pendulum experimental setup. (a)Mode1 vibration. (b)Mode2 vibration. . . . . . . . . . . . . . . . . . . . . . . . . . 117 5.2 NI Vision Builder on two degree of freedom experimental setup. . . . . . . . 119 viii 5.3 Experimental result on 2DOF inverted pendulum under periodic forcing am- plitude " = 0:0085 0:0005m and periodic forcing frequency ! = 30:0 0:1Hz. System is under mode 1 vibration. (a)periodic forcing (b)pendulum base horizontal movement (c)time traces of system response (d)time traces of system angular velocity (e)time traces zoom-in view of system response (f)time traces zoom-in view of system angular velocity. . . . . . . . . . . . . 122 5.4 Experimental result on 2DOF inverted pendulum under periodic forcing am- plitude " = 0:0085 0:0005m and periodic forcing frequency ! = 30:0 0:1Hz. System is under mode 1 vibration. Left graphs indicate the system response of the rst pendulum; right graphs indicate the system response of the second pendulum. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123 5.5 Power spectral density graphs of 2DOF inverted pendulum under periodic forcing amplitude " = 0:0085 0:0005m and periodic forcing frequency ! = 30:0 0:1Hz. System is under mode 1 vibration. (a)PSD for the rst pendulum (b)PSD for the second pendulum. . . . . . . . . . . . . . . . 124 5.6 Experimentalresulton2DOFinvertedpendulumunderperiodicverticalforc- ing. The forcing frequency was ! = 28:0 0:1Hz. System is under mode 2 vibration. Notice that the input forcing was not harmonic which had two forcing amplitudes involved. . . . . . . . . . . . . . . . . . . . . . . . . . . . 126 5.7 Experimentalresulton2DOFinvertedpendulumunderperiodicverticalforc- ing. The forcing frequency was ! = 28:00:1Hz. System is under mode 2 vibration. Left graphs indicate the system response of the rst pendulum; right graphs indicate the system response of the second pendulum. . . . . 127 5.8 Power spectral density graphs of 2DOF inverted pendulum under periodic vertical forcing. The forcing frequency was ! = 28:0 0:1Hz. System is under mode 2 vibration. (a)PSD for the rst pendulum (b)PSD for the second pendulum.. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 5.9 Experimentalresulton2DOFinvertedpendulumunderperiodicverticalforc- ing with 0 initial condition. The forcing frequency was ! = 30:0 0:1Hz. Notice that the input forcing was not harmonic which had two forcing am- plitudes involved. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 129 5.10 Experimentalresulton2DOFinvertedpendulumunderperiodicverticalforc- ing with 0 initial condition. The forcing frequency was ! = 30:0 0:1Hz. Notice that the input forcing was not harmonic which had two forcing am- plitudes involved. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 131 ix 5.11 Power spectral density graphs of 2DOF inverted pendulum under periodic vertical forcing with 0 initial condition. The forcing frequency was ! = 30:00:1Hz. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132 6.1 Modi ed 2DOF model has each pendulums center of mass o¤set to a new location di¤erent than the pivot point. . . . . . . . . . . . . . . . . . . . . . 134 x Abstract A pendulum is statically unstable in its upright inverted state due to the Earths gravitionalattractionwhichpointsdownward. However, withproperforcing, thependulum can be stabilized in its upright inverted state. Special interest is on periodic vertical forc- ing applied to the pendulums base to stabilize it around the upright inverted equilibrium. Many researchers have studied how to stabilize the system by varying various parameters, in particular its amplitude and frequency. Most have focused on the single degree of free- dom inverted pendulum case, which with linear assumption can be described via Mathieus equation. The system stability can then be characterized by Floquet theory. Our focus is on searching for the periodic solutions inside the linearly stable region of the pendulums invertedstatewhenthependulumisunderproperperiodicforcing. Ourresearchshowsthat under appropriate excitation by controlling the forcing amplitude and frequency, the pen- dulum can maintain certain periodic orbits around its inverted state which we characterize in a systematic way. Inthisthesis,weappliedfourdi¤erentkindsofgeometricrealizationsofthesystem response: systemtimetraces,systemphaseportraits,threedimensionalviewsofthesystem phase portrait as a function of input forcing, and the systems power spectral density dia- gram. By analyzing these four diagrams simultaneously, we characterize di¤erent kinds of multi-frequencyperiodicbehavioraroundthependulumsinvertedstate. Tofurtherdiscuss the e¤ect of the nonlinearity, we applied perturbation techniques using the normalized forc- ing amplitude as a perturbation parameter to carry out the approximate periodic solutions on a single degree of freedom inverted pendulum nonlinear model. xi We also discuss the multiple degree of freedom inverted pendulum system. Both numerical simulation and experiments were performed and detailed comparisons are dis- cussed. Our numerical simulations show close qualitative agreement with experiments. 1 Chapter 1 Inverted Pendulum Dynamics 1.1 Introduction This dissertation discusses the e¤ect of pendulum stabilization around its inverted state under vertical high-frequency periodic forcing. Although, in general, the inverted pendulums response wont be periodic under periodic forcing, we are interested in charac- terizing periodic and multi-periodic behavior. Of special interest is the inverted pendulum stabilizedinperiodicorbitswherethesystemperformsmultiplefrequencynoddingbehavior as discussed in Acheson[1], [2]. Thisphenomenon, is alsoknown as hoveringmotion(refer toWeibelandBaillieul[60]). Thesystemsphaseportraitexhibitsperiodicorbitswhichare similar to limit cycles to a time-invariant system. In this dissertation, we study the system through four di¤erent types of diagrams: system time traces, system phase portraits, three dimensional views of the system phase portrait as a function of input forcing, and the sys- tems power spectral density diagram. First, we discuss the single degree of freedom system showing modeling, numerical simulation and experiments. In later chapters, we extend the 2 single degree of freedom system to higher degree of freedom cases. Our system response resultsshowcloseagreement(qualitatively)betweennumericalsimulationandexperiments. Beginningfromthesingledegreeoffreedominvertedpendulummodelunderverti- cal periodic forcing, a set of nonlinear non-autonomous ordinary di¤erential equations were obtained to describe our model. First, we follow those remarkable ancestors in researching the stability in the single degree of freedom inverted pendulum problem such as Acheson [1], [2], Weibel and Baillieul [60], Rand and Morrison [46], [38]. The linearized equation of motion (Mathieus equation) is shown in early sections of this dissertation. We apply Floquet theory to Mathieus equation, which is also known as Hills equation, to obtain the stability diagram through transition curves calculated from Hills equation. To under- stand the e¤ect of the nonlinearity, we then apply a two-timing perturbation technique by assigning a (normalized) forcing amplitude as the perturbation parameter. We then carry out analytic solutions from the nonlinear single degree of freedom model. Later, we perform numerical simulations based on the nonlinear single degree of freedom model to obtain periodic solution sets which we characterize using di¤erent normalized forcing para- meters and initial conditions. Finally, we conducted a sequence of experiments to obtain detailed measurements to compare with our simulation results. Although experimental re- sults show damping e¤ects, which contradict our assumption of energy conservation, for short times and for the purposes of identifying the multi-frequency state, the experiments and simulations matched reasonably well. Presented in this paper is a way to analyze a dynamic system where the non- linearity and time-varying properties cannot be neglected. Unlike the well-known simple 3 pendulum example in classical mechanics, the inverted pendulum doesnt provide the same kind of clockwork regularity which can be properly predicted by the solutions of a simple timeinvariantlinearordinarydi¤erentialequation. Furthermore,suchregularityinthesim- ple pendulum can be properly modelled by conservation laws where the system is regular back-and-forth swinging motion can be seen from exchanging its potential and kinetic ener- gies. On the other hand, the inverted pendulum needs certain external energy to excite the systemtoovercometheEarthsgravitationalattractioninordertostabilizeitinastatically upright state. The periodic excitation to a simple pendulum brings additional complexity to the system and its governing equation of motion usually needs to be treated as non- autonomous, non-conservative and nonlinear. The analysis of such systems is a challenging task because there is no closed-form analytical solution that can be directly applied. We use a 4th order Runge-Kutta numerical integration technique to get solutions fordi¤erentinitialconditions. Byexaminingthetimeevolutioninitsphaseportrait,wecan nd limit cycles under speci c driving parameters and initial conditions. In later chapters, we apply singular perturbation techniques to the nonlinear model to derive sequences of asymptotic approximations to get higher order of accuracy to the solution of the inverted pendulum under multi-frequency nodding state which can be compared directly to the numericalintegrationandexperimentalresults. Perturbationanalysisalsoprovidesabetter understanding of the complex system and stimulated the development of new methods for the numerical solution of the higher order approximations. Experiments were then introduced and compared with the analytical and numerical integration results in detail. 4 1.2 Motivation and Related Early Work on Inverted Pendu- lum Dynamics The instability of the pendulum in its inverted state can be removed by applying rapid vertical oscillatory forcing to its pivot point. This phenomenon had been discovered in early 20th century by Stephenson [53] [54] . When placing it initially in its upright state, the inverted pendulum will remain in its upright vertical state without falling over. However, this is not the only stable state when pendulum is inverted. Acheson et al [1] [2] [3] discovered that the pendulum may withstand large initial disturbances when applying appropriateverticaloscillatoryforcingwithproperdrivingamplitudeandfrequency. Under such conditions, the pendulum will oscillate around its upright vertical equilibrium point without falling over. The pendulum becomes trapped in a limit cycle oscillation. This is an interesting phenomenon which can be applied to model many engineering systems undergo- ingperiodicvibrations. Manyresearchpapersdiscusstheoscillationsofadrivenpendulum. Here are some examples: Blackburn et al [12] studied the stability and Hopf bifurcations in an inverted pendulum; Smith et al [49] investigated the behavior of an inverted pendulum through experimental measurements; Kalmus [25] worked on a driven inverted pendulum experimentusingaspeaker;Michaelis[36]usedanelectricjigsawtodriveaninvertedpendu- lum and study its behavior through stroboscopic photos; Acheson et al [1],[2],[3] compared thestabilityofaninvertedpendulumfromtheoreticalmodelsandexperimentalapproaches; Flashner et al [18] studied the bifurcation through the point mapping method. This kind of system has a time-dependent coe¢ cient in a governing di¤erential equation also known as parametric excitation. The physical sense of the parametric excitation is that the sti¤ness 5 ormomentofinertiainadynamicsystemdependsontime. Thisisquitedi¤erentcompared with adding a nonhomogeneous term in the governing di¤erential equation which results from external excitation. Parametric excitation occurs in a wide variety of engineering applications. Here are a few recent applications: Stephan G. et al [51] investigated para- metricexcitationinhigh-speedmillingapplications; Gani[19]studiedparametricexcitation of stay-cables; Yu realized parametric excitation in a nanowire system using an oscillating electric eld; Mennem [35] studied parametrically excited vibrations in spiral bevel geared systems; Kaajakari [26] realized parametric excitation in ultrasonic surface micromachine actuation. Mostanalysismethodsfortime-varyingnonlineardynamicsystemsarecoordinate- based approaches. This dissertation aims to develop a systematic technique to understand andanalyzetime-varyingand, morespeci cally, time-periodicnonlinearinvertedpendulum dynamic systems in a geometric setting. Series expansion methods and averaging theory were two powerful tools used in solving time varying nonlinear di¤erential equations. Vela [64],[65] gives a detailed treatment of the method of averaging and the related theorems that comprise averaging theory which showed how averaging theory may be used to time- dependent di¤erential equations and it can also be applied to nonlinear control theory of underactuated systems. He pointed out that the averaging theory is derived from nonlinear Floquettheoryandperturbationtheory. Velasworkalsoshownthatnonlineartime-varying control systems can also be given an exponential representation, meaning that all of the intuition and analysis from linear control theory may provide the control engineer with the needed background to construct and analyze stabilizing controllers for nonlinear systems. 6 Bogoliubov and Mitropolsky [13] cover the same topics, and also give a general algorithm for calculating higher-order averages. The Poincaré map is a powerful tool for visualizing the results from Floquet theory and extending it to the nonlinear regime. Bogoliubov and Mitropolsky were aware that the averaged equations of a time-dependent di¤erential equation gave the Poincaré map, thereby allowing for stability analysis, and also posited that their method was able to recover higher orders of averaging. This method of averaging is known as the Krylov- Bogoliubov-Mitropolsky (KBM) method of averaging. The higher- ordermethodsproposedbyBogoliubovandMitropolskyhavebeenstudiedandextendedby several researchers as they are the most general and most powerful of the known averaging methods. Averaging theory also includes two-timing methods, which involve fast and slow time scale; only the fast time scale is averaged over. In many cases averaging over multiple dimensions will introduce resonance. There also exist other papers detailing higher-order averaging theory, however they focus on particular classes of time-periodic systems. 7 Chapter 2 Inverted Pendulum Under Periodic Vertical Forcing 2.1 Single Degree of Freedom Inverted Pendulum Model Equations of Motion Considerasimpli edsingledegreeoffreedominvertedpendulumwhichconsistsof a weightless rigid rod with lengthL. A point massm is mounted on one end of the inverted pendulum and on the other end a frictionless pivot point P is subjected to periodic vertical forcing with coordinate u(t) ="cos!t measured downward from some xed point O. Here " denotes the vertical driving amplitude and ! denotes the periodic driving frequency. The whole system is under uniform gravitational eld with gravitational acceleration g pointing downward. The pendulum has angular displacement measured from the vertical upright position(inverted state) as shown in Figure 2:1. Here we neglect the e¤ect of damping in 8 order to apply the energy conservation property by means of Lagranges equations. The equations of motion can be derived from the general form of Lagranges equations as d dt @T @ _ @T @ + @V @ = (2.1) d dt @T @_ u @T @u + @V @u = U (2.2) where T is the kinetic energy, V is the potential energy and and U are generalized nonconservative forces with respect to and u direction respectively. Since the vertical drivingvibrationistheonlysourceofnonconservativeforces,therighthandsideofequations (2:1) and (2:2) can be written as = 0 , U =F: (2.3) From the energy conservation law, (2:1) indicates the energy exchange between kinetic energy and system potential due to the rotational motion of the pendulum; (2:2) indicatesthebalanceofverticalperiodicforcingandsystemenergytransferbetweenkinetic energy and potential energy due to the applied force. The kinetic energy can be derived as T = 1 2 m v 2 x +v 2 y = 1 2 m L _ cos 2 + _ u+L _ sin 2 = 1 2 m L 2 _ 2 +2L_ u _ sin+ _ u 2 (2.4) where v x is horizontal velocity vector and v y is vertical velocity vector of the point mass as shown in Figure 2:1. The potential energy V depends on the vertical displacement y 9 Figure 2.1: Single degree of freedom inverted pendulum under vertical periodic forcing at its base P: which can be written as V =mgy =mg[L(1+cos)+u]: (2.5) As indicated in (2:5), the system potential only depends on the posture of the pendulum at that instance, not its motion. In other words, V is not a function of _ or _ u. We can then calculate individual components of (2:1) from derivatives of (2:4) and (2:5) as @T @ _ = mL 2 _ +mL_ usin (2.6) d dt @T @ _ = mL 2 +mL usin+ _ u _ cos @T @ = mL_ u _ cos @V @ = mgLsin 10 From derivatives of (2:4) and (2:5) to individual components of (2:2) we get @T @_ u = m L _ sin+ _ u (2.7) d dt @T @_ u = m L sin+L _ 2 cos+ u @T @u = 0 @V @u = mg Hence, substitute (2:6) and (2:7) into (2:1) and (2:2), we obtain the explicit Lagranges equations as mL 2 +mL usinmgLsin = 0; (2.8) m L sin+L _ 2 cos+ u +mg =F: (2.9) From time derivatives to the driving displacement u , we obtain its velocity as _ u = "!sin!t. From second time derivatives to u, we obtain its acceleration as u ="! 2 cos!t: (2.10) The simpli ed equation of motion of the single degree of freedom inverted pendulum model under vertical vibration can then be obtained by inserting (2:10) into (2:8) as mL 2 mgL+mL"! 2 cos!t sin = 0: (2.11) 11 Equation (2:11) is a nonlinear di¤erential equation with a time-dependent coe¢ cient. By dividing both sides of (2:11) by mL 2 we obtain g L + "! 2 L cos!t sin = 0: (2.12) Unlikethesimplependulumcasewhosestatecoe¢ cientisconstant,oursystemhas a time-varying coe¢ cient due to driving forces. When the driving acceleration amplitude is smallerthangravity,"! 2 <g,thesystemcoe¢ cient, g L + "! 2 L cos!t;doesntchangesignand thesystemresponsehaslittlee¤ectduetoforcing. Basedonlineartheory,thestateissaidto behyperbolicsincetheeigenvaluesoftheJacobianmatrixaroundtheequilibriumpointsdo notlieontheimaginaryaxis. Inthiscase, thenonlinearbehaviorcanbeproperlypredicted through linear analysis around its equilibrium points. On the other hand, when the driving acceleration amplitude is larger than gravity, "! 2 > g, the system coe¢ cient changes sign. In this case, the equilibrium point will change from hyperbolic to non-hyperbolic since the eigenvalues of the Jacobian matrix around its equilibrium points will sometimes lie on the imaginary axis. De nition 1 For a continuous-time system, an equilibrium is called non-hyperbolic if the Jacobian matrix evaluated at the equilibrium point has at least one eigenvalue with zero real part. In other words, if there is any eigenvalue located on the imaginary axis then the system is non-hyperbolic. For a discrete-time system a xed point is called non-hyperbolic if the Jacobian matrix evaluated at the xed point has at least one multipliers located on the unit circle. The nonlinear behavior of a hyperbolic system can be properly predicted from its 12 linearized counterparts. The Hartman-Grobman theorem states that the local phase por- trait near a hyperbolic equilibrium point of the nonlinear system is topologically equivalent to the phase portrait of its linearization. One important consequence is the stability type of the equilibrium point of the nonlinear system is preserved by its linearization. However,if the system is non-hyperbolic, the theorem does not hold. In this paper, we are exclusively discussingthee¤ectofdrivingaccelerationamplitudemuchlargerthangravityor"! 2 >>g which will excite the pendulum to stabilize at its inverted state and such a system usually is non-hyperbolic. Acheson [3],[1] shows that the pendulum can be stabilized in its inverted state only when " 2 ! 2 > 2gL which guarantees "! 2 >> g under small driving amplitude ". Our experimental results also show that a driving acceleration of about 30G 50G is required for a single degree of freedom pendulum to stabilize its inverted state, i.e. "! 2 needs to be 30 to 50 times larger than g. The analysis of time-varying, nonautonomous nonlinear systems is a challenging task. The time-varying coe¢ cient or parametric excitation of (2:12) plays a big role in contributing to the response of the system when the driving acceleration is much larger than the gravitational acceleration, or "! 2 >> g. When the frequency of excitation is su¢ ciently far from the primary resonance, a small parametric excitation can produce a large response. This time-varying property makes the system di¢ cult to solve analytically and the system response is highly dependent on the initial conditions. However, the system response can be simpli ed if we limit the time-varying parameters to be time-periodic. Withtime-periodicparametricexcitation,theremayexistperiodicsolutions(fromaveraging theory)ifthesysteminitialconditionsareselectedproperly. Averagingtheoryisatechnique 13 which synthesizes Floquet theory and perturbation theory by applying series expansions to approximate time-varying vector elds [64], [65]. In the next section, we apply averaging theory to a linear single degree of freedom inverted pendulum system by means of Floquet theory and series expansion through system perturbation parameter. System stability can be shown using transition curves and Stutts diagram. Linearized sDOF Model and Mathieu s Equation When using the small angle assumption, sin ! , the original nonlinear single degree of freedom inverted pendulum equation of motion (2:12) can be linearized as g L + "! 2 L cos!t = 0: (2.13) Equation (2:13) is then transformed to well-known Mathieus equation d 2 d 2 +(+cos) = 0: (2.14) Here, the rst system parameter = !n ! 2 is de ned as negative value of the square of the ratio of the system natural frequency ! n = q g L and the driving frequency !. The secondsystemparameter = " L isde nedasnormalizeddrivingamplitude. Equation(2:14) has new system time de ned as =!t. Primary Resonance of the Transition Curves of Mathieus Equation Apply the 2-timing perturbation technique by introducing two new variables: fast time = and slow time = to (2:14). The systems state is then translated into ()7!(;) and its derivatives can be carried out using chain rule as follows: 14 _ = d d = @ @ d d + @ @ d d = @ @ + @ @ ; (2.15) = d 2 d 2 = @ 2 @ 2 +2 @ 2 @@ + 2 @ 2 @ 2 : (2.16) Substitute (2:16) into (2:14), @ 2 @ 2 +2 @ 2 @@ + 2 @ 2 @ 2 +(+cos) = 0: (2.17) An approximate solution of (2:14) can be obtained from expanding in a power series for small as: (;) = 0 (;)+ 1 (;)+ 2 2 (;)+::: (2.18) Substituting (2:18) into (2:17) and neglecting higher order terms of O 2 we get @ 2 ( 0 + 1 ) @ 2 +2 @ 2 ( 0 + 1 ) @@ + 2 @ 2 ( 0 + 1 ) @ 2 +(+cos)( 0 + 1 ) = 0 (2.19) Collect terms of the same orders of : @ 2 0 @ 2 + 0 = 0; (2.20) @ 2 1 @ 2 + 1 = 2 @ 2 0 @@ 0 cos: (2.21) The general solution to (2:20) can be expressed as 0 (;) =A()cos p +B()sin p : (2.22) Di¤erentiating to (2:22) with respect to and we have @ 2 0 @@ = p dA() d sin p + p dB() d cos p : (2.23) 15 Substituting (2:22) and (2:23) into (2:21) , we then obtain @ 2 1 @ 2 + 1 = 2 p dA() d sin p + p dB() d cos p A()cos p +B()sin p cos: (2.24) Apply the following trigonometric identities to simplify (2:24) sincos = 1 2 (sin(+)sin()); (2.25) coscos = 1 2 (cos(+)+cos()); (2.26) where , are dummy variables. Substitute (2:25) and (2:26) using proper variables in (2:24) to replace and we obtain @ 2 1 @ 2 + 1 = 2 p dA() d sin p 2 p dB() d cos p A() 2 cos p +1 +cos p 1 B() 2 sin p +1 sin p 1 : (2.27) For general values of , removal of resonance terms require the coe¢ cients of sin p and cos p become zero which gives the trivial solutions as: dA() d = 0; (2.28) dB() d = 0: (2.29) This means that for general ; the cos!t driving term in Mathieus equation (2:14) has no e¤ect. However, if we choose = 1 4 and substitute into (2:27) we get 16 @ 2 1 @ 2 + 1 4 1 = dA() d sin 2 dB() d cos 2 A() 2 cos 3 2 +cos 2 B() 2 sin 3 2 sin 2 : (2.30) From analyzing the right-hand side of (2:30), the removal of resonance terms gives the following solutions: dA() d = B() 2 ; (2.31) dB() d = A() 2 : (2.32) We may also rewrite (2:31) and (2:32) in matrix form as 2 6 6 4 dA() d dB() d 3 7 7 5 = 2 6 6 4 0 1 2 1 2 0 3 7 7 5 2 6 6 4 A() B() 3 7 7 5 : (2.33) Combine (2:31) and (2:32) to obtain d 2 A() d 2 = A() 4 : (2.34) Here A() and B() involve exponential growth where the instability occurs when = 1 4 : This corresponds to a 2 : 1 subharmonic resonance in which the driving frequency is twice the natural frequency [38]. By expanding in a power series in we may get a curve due to resonance excitation in more generalized way as = 1 4 + 1 + 2 2 +::: (2.35) Substitute (2:35) into (2:19) and neglecting terms of O 2 we get 17 @ 2 ( 0 + 1 ) @ 2 +2" @ 2 ( 0 + 1 ) @@ +" 2 @ 2 ( 0 + 1 ) @ 2 + 1 4 + 1 +cos ( 0 + 1 ) = 0: (2.36) Collect terms of the same orders, yields @ 2 0 @ 2 + 1 4 0 = 0; (2.37) @ 2 1 @ 2 + 1 4 1 = 2 @ 2 0 @@ 0 cos 1 0 ; (2.38) which results in the following additional terms in (2:31) and (2:32) as dA() d = 1 1 2 B(); (2.39) dB() d = 1 + 1 2 A(): (2.40) We may also rewrite (2:39) and (2:40) in matrix form as 2 6 6 4 dA() d dB() d 3 7 7 5 = 2 6 6 4 0 1 1 2 1 + 1 2 0 3 7 7 5 2 6 6 4 A() B() 3 7 7 5 ; (2.41) Comparing with (2:34); an additional term is introduced d 2 A() d 2 + 2 1 1 4 A() = 0: (2.42) From (2:42) we may conclude that if 2 1 1 4 > 0 , that is, if either 1 > 1 2 or 1 < 1 2 , then A() and B() will be sine and cosine function of slow time . We may substitute boundary values of 1 into (2:35); which yields = 1 4 2 +O 2 : (2.43) 18 Therearetwocurvesdescribedin (2:43);startingfromthepoint = 1 4 onthe axisin plane. These are recognized as primary resonances of the transition curves [46]. Transition curves as in Figure 2.2 represent stability changes, region in between of these two curves is unstable. Inside the unstable region, for small , grows exponentially in time and the solution is unbounded. Outside the unstable region, from (2:22), (2:39), (2:40) and (2:42), is the sum of terms each of which is the product of two periodic or harmonic functions and is bounded, speci cally, (t) is a quasiperiodic function of time. Floquet Theory Because the time-varying coe¢ cient in Mathieus equation (2:14) is periodic in time, we may apply Floquet theory to simplify the system. Floquet theory is concerned with the following system of rst order di¤erential equations: dx dt =A(t)x; (2.44) where x is an n1 column vector, and A is an nn time-periodic coe¢ cients matrix with period T with the following property: A(t+T) =A(t): (2.45) Although A(t) varies periodically in time, the general solutions of (2:44) are typically not periodic. Floquet theory indicates that there exists a fundamental matrix solution A 0;t of (2:44) of the form A 0;t =P (t)exp(t); (2.46) 19 where P (t) is period in time P (t+T) =P (t), and is in general a complex number. The general solution of (2:44) takes the following form: x(t) = n X i C i exp( i t)P i (t) (2.47) where C i are constants that depend only on initial conditions x(t 0 ) = x 0 : They take the following form: C = A 0;t 0 1 x 0 : P i (t) in (2:47) are vector-valued functions which are also periodic in time with period T, and i are complex numbers known as Floquet exponents. The Floquet multipliers can be de ned as i = exp( i T): (2.48) The Floquet exponents are not unique since exp i + 2ik T T = exp( i T) for any integer value k . Furthermore, the long-term behavior of the solution of x(t) is determined by the Floquet exponents. Floquet theory allows us to reach an important conclusion about the solutions stability. If any of the Floquet multipliers have modulus greater than unity, i.e. j i j > 1, then x(t) is unbounded as t ! 1 and the system becomes unstable. On the other hand, if all Floquet multipliers have modulus less than unity, i.e. j i j < 1, then x(t) is bounded as t ! 1 and the system is stable. The Floquet exponents and Floquet multipliers represent the growth rate of di¤erent perturbations averaged over a cycle. Floquet exponents are rates with unit t 1 and Floquet multipliers are dimensionless numbers that give the period to period increase or decrease of the perturbation. 20 Floquet theory can be applied to a generalized Mathieu equation (also known as Hills equation [46]) which takes the following form: d 2 x dt 2 +f(t)x = 0; f(t+T) = f(t): (2.49) Here x and f are scalars, and f(t) represents a general periodic function with period T . Thesecondorderordinarydi¤erentialequationin (2:49)canbetransferredintoasystemof two rst order o.d.es by de ning x 1 =x and x 2 = _ x and then substitute into (2:49) yields 2 6 6 4 _ x 1 _ x 2 3 7 7 5 = 2 6 6 4 0 1 f(t) 0 3 7 7 5 2 6 6 4 x 1 x 2 3 7 7 5 : (2.50) Two fundamental solutions can then be constructed as 2 6 6 4 x 11 (t) x 12 (t) 3 7 7 5 and 2 6 6 4 x 21 (t) x 22 (t) 3 7 7 5 which satisfy the initial conditions 2 6 6 4 x 11 (0) x 12 (0) 3 7 7 5 = 2 6 6 4 1 0 3 7 7 5 and 2 6 6 4 x 21 (0) x 22 (0) 3 7 7 5 = 2 6 6 4 0 1 3 7 7 5 respectively. We can then construct a matrix C which is the fundamental solution matrix evaluated at time T C = 2 6 6 4 x 11 (T) x 12 (T) x 21 (T) x 22 (T) 3 7 7 5 : (2.51) From Floquet theory, we may conclude that the system stability is determined by the eigenvalues of C. If all eigenvalues of C have modulus less than unity, the systems solution will be bounded and the system is stable; if any eigenvalue of C is greater than unity, the systems solution will grow exponential in time and the system is unstable. The eigenvalues of C can be written as 2 tr(C)+det(C) = 0; (2.52) 21 where tr(C) denotes the trace of C which de nes as tr(C) =x 11 (T)+x 22 (T) and det(C) denotesthedeterminantofC whichde nesas det(C) =x 11 (T)x 22 (T)x 12 (T)x 21 (T). The Hills equation (2:49) has a special property that its determinant has value equal to one, i.e. det(C) = 1. Due to this special property,(2:52) becomes 2 tr(C)+1 = 0; (2.53) and the eigenvalue has solution = tr(C) q (tr(C)) 2 4 2 : (2.54) Therefore, theeigenvaluehasrealrootswhenjtr(C)j> 2andtheeigenvalueshave apairofcomplexconjugaterootswhenjtr(C)j< 2. Anotherspecialpropertyistheproduct of the eigenvalues must equal to unity. Thus, in the case of two real roots, if one root is less than unity then the other root must be larger than unity which results in instability. There are two stable conditions for two real roots case which can be obtained whenjtr(C)j = 2: If tr(C) = 2 then two real eigenvalues are 1;2 = 1;1 and according to Floquet theory, this condition corresponds to a periodic solution with period T: If tr(C) = 2 then two real eigenvalues are 1;2 = 1;1 and according to Floquet theory this condition corresponds to a periodic solution with period 2T: On the other hand, in the case of a pair of complex conjugateroots,botheigenvaluesmustlieontheunitcircle. Wecanconcludethiscondition asneutralstabilityandthesystemhasquasiperiodicsolution. AllstableconditionsinHills equation have eigenvalues which lie on the unit circle. Such a system is non-hyperbolic. Since Hills equation has no damping e¤ect and its stable conditions are non- hyperbolic, it allows us to make conclusions about long time behavior using the solution 22 formmerelyoneforcingperiod. Fulltransitioncurvescanbecarriedoutusingtheharmonic balance technique by applying the fact that the period of the forcing function in Mathieus equation (2:14) is 2 periodic, i.e. T = 2: Detail discussion about Hills equation and harmonic balance is shown in Rand [46]. Transition Curves of Mathieus Equation: Strutts Diagram The transition curves on the plane also known as Strutts diagram is shown in Figure 2.2. The transition curves were obtained by applying Floquet theory to linear Mathieus equation and then applying the harmonic balance technique. Refer to the Ap- pendix at the end of this paper for a listing of the complete sets of transition curves derived by Rand [46]. In Figure 2.2, the region of > 0 denotes the pendulum has downward vertical equilibrium where inside the bounded region(shaded area) the pendulum will be trapped arounditsdownwardverticalposition(normalstate); theregionof < 0denotesthependu- lumhasuprightverticalequilibriumwhereinsidetheboundedregion(shadedarea)thepen- dulum will be trapped around its upright vertical position(inverted state). The unbounded regionbetweenthetwotransitioncurveswereobtainedthroughresonanceexcitation, start- ing from = n 2 4 ;n = 0;1;2;3::: (2.55) for small . The unbounded region from the starting point = 1 4 at n = 1 is called the primary resonance and its transition curves were calculated from (2:43). The unbounded region under primary resonance has the largest area. As n value increases, the unbounded region moves to the right hand side of the primary resonance and its covering area becomes 23 Figure 2.2: Transition curves for the linear Mathieu equation on the plane separate regions of bounded solution(shaded area) from regions of unbounded solution or unstable system(white area). 24 smaller. We only discuss the region of < 0 when the pendulum has upright vertical equilibrium. In other words, our special interest focuses on the bounded region(shaded area) of < 0 covered in between the two transition curves of n = 0 and n = 1. Figure 2.3 shows the zoom-in view of Figure 2.2 where < 0 region being emphasized. Recall that in our inverted pendulum model (2:14) was de ned as negative value of the square of the ratio of the system natural frequency and the driving frequency which takes the following equation: = !n ! 2 . If systems natural frequency remain xed then the higher the forcing frequency !, the closer approaches 0. There are 3 regions labelled A;B;C in Figure 2.3 which indicates our region of interest in our single degree of freedom inverted pendulum model. These regions are in relative higher forcing frequency region (> 25Hz) where is close to 0. In our simulation trials, we xed the forcing frequency while varying the forcing amplitude gradually from A to B then to C regions. Our main focus is on B region where the pendulum is stabilized in its inverted state, while A and C regions are only to test stability boundaries. Another way to present the stability boundaries on a graph of an inverted pendu- lum under periodic forcing is to graph the stability region in the ! plane as in Figure 2.4. The shaded area in Figure 2.4 indicates the results from Acheson and Mullin [3] which showed the stability boundaries calculated from solving linear Mathieus equation. Inside the shaded region, the pendulum can be stabilized around its inverted posture where the system has bounded solutions. On the other hand, the region outside the shaded region has unbounded solutions. Left(n = 0, blue-color curve) and right(n = 1, green-color curve) transitioncurvesin ! stabilitydiagraminFigure2.4showthestabilityboundariesofour 25 Figure 2.3: Stutts diagram for stability regions separated by the transition curves of n = 0 and n = 1 within the range of0:5< < 0:5. 26 calculations matched with the results from Acheson and Mullin [3]. These two transition curves are the same curves as shown in Figure 2.3 inside < 0 region. The left stabil- ity limit curve described by " 2 ! 2 = 2gL was calculated from Acheson [1] which showed more conservative stability region compared with the left transition curve(n = 0). This left stability curve(red-color) is very close to the left transition curve(blue-color) as shown in Figure 2.4. Those green-color data points denote the experimental results from Acheson and Mullin [3] which show the stability boundary from their experimentations were even wider than their prediction from solving linear sDOF inverted pendulum model. This can be explained by the e¤ect of system damping due to friction in their experimental setups. In other words, system damping helps the systems stability in this case. There are three sets of diamond data points in Figure 2.4, indicating 3 sets of periodic solutions from our numerical simulations to the nonlinear single degree of free- dom inverted pendulum model. Basically, these diamond data points indicate the sta- bility boundaries and periodic solutions to the same single degree of freedom inverted pendulum model where we xed the normalized forcing frequency ! at three locations: ! = 11:6218; ! = 14:2784; and ! = 20:0, while we varied normalized forcing amplitude from small to large values. Our simulation results showed that the left limit was located on the left transition curve, while the right limit was wider than the right transition curve. The periodic solutions were discovered in numerical simulation inside the stable region in Figure 2.4 labelled as 4 : 1; 5 : 1; 6 : 1;..., where the ratio indicates the forcing period with respect to systems periodic response. We denoted this ratio as frequency ratio N r whose detail will be shown in later chapters. The set of data points in blue-color circles 27 Figure 2.4: Stability regime diagram of single degree of freedom inverted pendulum under vertical periodic forcing showing in the ! plane. 28 indicates the parameter distributions of our single degree of freedom inverted pendulum ex- perimentations. Due to physical limitations of our experimental setups, experimental data only occupied the region of normalized forcing amplitude < 0:2 and normalized forcing frequency 10 ! 20. More details of Figure 2.4 will be shown in later sections in this paper. Numerical Integration of the Nonlinear sDOF Inverted Pendulum Model Simulation Based on Parameter Continuation Code Using 4 th Order Runge- Kutta Numerical Integration Scheme This section discusses the simulations from direct numerical integration of the time-varying single degree of freedom inverted pendulum model. We normalized the sDOF modelbysubstitutingthenormalizeddrivingamplitudeandnormalizeddrivingfrequency ! into (2:12). The normalized driving amplitude is de ned as = " L and normalized driving frequency is de ned as ! = ! !n , where the system natural frequency is ! n = q g L . Therefore, (2:12) can be expressed as: (t) 1+ ! 2 cos!t ! 2 n sin(t) = 0: (2.56) The original 2 nd order system equation can then be re-written as two 1 st order equations _ x 1 = x 2 ; (2.57) _ x 2 = (1+ ! 2 cos!t)! 2 n sinx 1 : where x 1 = ; and x 2 = _ : The natural frequency of the system is ! n = q g L ; the normalized driving amplitude is = " L , and the normalized driving frequency is ! = ! !n . 29 The periodic forcing in the system can be replaced by introducing two more dimensions as _ u = u(1u 2 v 2 )!v; (2.58) _ v = v(1u 2 v 2 )+!u: where the two dynamic harmonic functions are described as u = cos!t and v = sin!t. We apply a 4 th order Runge-Kutta numerical integration scheme to (2:57) and (2:58). Four Geometric Representations of the System Response In this paper, we apply four di¤erent kinds of geometric realizations of system response: system time traces, system phase portrait, three dimensional view of system phase portrait compare to input forcing, and systems power spectral density diagram. Below are short descriptions to these four representations. System Time Trace The system time trace represents the system response with respect totimeevolution. Periodicresponsesofasystemcanbedetectedbyinspectingtherepeated patternsinasystemtimetrace. Ifasystemresponserepeatsitselfaftersomeconstanttime T, or x(t+T) =x(t), then such a system is periodic with period T. SystemPhasePortrait Aphaseportraitisageometricrepresentationofthetrajectories of a dynamic system in the phase plane. A system with periodic response exhibits closed trajectories or orbits in a phase portrait. Associated with our phase portrait graphs are vector elds showing piecewise system response in a phase plane starting from uniformly distributed initial conditions. 30 3DViewofSystemPhasePortrait Byaddinganewdimension,thesystemsperiodic input u = "cos!t; to a system phase portrait, we may get a three dimensional view of a dynamicresponse. Allclosedorbitsshowninaperiodicphaseportraitgraphresultinclosed curves in a 3D phase portrait. Recall that the system response in a periodically forced inverted pendulum system is, in general, non-periodic. However, if a systems response is periodic under periodic forcing, one period of the system response cycle must be an integer multiple of the periodic forcing period. In other words, the periodic forcing frequency is N r times larger than the periodic system responses frequency where N r is an integer. Therefore, the frequency ratio N r can be measured from counting the cycles of a 3D phase portrait graph of a periodic response. Power Spectral Density (PSD) Analysis Power spectral density(PSD) describes how the power of a signal or discrete time series is distributed with frequency. PSD requires that the Fourier transforms of the signals exist and that the signals are square-integrable or square-summable. PSDcanbeappliedtoextractthefrequencycomponentsoftheresponse in order to obtain the periodic pattern of the signals composed of harmonic functions. We apply the discrete Fourier transform(DFT) to the dynamic response as follows: u k u(t k )u(kt) =a 0 +2 N=2 X n=1 a n cos 2nt k T +b n sin 2nt k T : (2.59) where k = 1;2;:::;N is an integer index of t k kt and t is the sampling interval; N is the total number of samples and T =Nt is the sample period. The sampling frequency is de nedasf s = 1 t whichleadstothemaximumobservablefrequency, orNyquistfrequency, 31 f ny = f s =2 or f ny = 1 2t . For convenience, throughout this paper for PSD calculations we used a uniform sampling frequency which is 500Hz. The spectral coe¢ cients a n and b n are de ned as a n = 1 N N X j=1 u j cos 2nj N ; (2.60) b n = 1 N N X j=1 u j sin 2nj N : (2.61) The spectral coe¢ cients represent the amplitudes of the harmonic components extracted from the dynamic response. The power spectral density is then de ned as: PSD =U U =a 2 n +b 2 n ; (2.62) where U is the spectrum of u(t) de ned as U n T a n ib n = 1 N N X j=1 u j exp i 2nj N : (2.63) where i p 1 and U represents the complex conjugate of U. Numerical Simulation Results of a sDOF Inverted Pendulum Webeginbypresentingthe rstthreetypesofsystemrepresentations: timetraces of inverted pendulum angles and angular velocities, phase portraits of the system, and a 3-dimensional view of the system phase portrait with respect to periodic forcing. Later in this chapter we introduced the power spectral density to further inspect its harmonic components. 32 Figure 2.5: Numerical integration results of the single degree of freedom inverted pendulum undernormalizedforcingamplitude = 0:1708andnormalizedforcingfrequency ! = 8:7313 with initial condition 0 = 0:1 and _ 0 = 0. (a)System time response. (b)Angular velocity time trace. (c)System phase portrait. (d)3D view of phase portrait with respect to periodic forcing. 33 From inspection of Figure 2.5(a) we see the system time trace repeats itself in a period around 1:5s; which refers to slow time of the system response. The fast time response is proportional to the periodic vertical forcing at normalized forcing parameters = 0:1708 and ! = 8:7313. From Figure 2.5(b) we can see that the system angular velocity response also has the same slow time period of around 1:5s: From the phase portrait graph in Figure 2.5(c), we see a symmetric periodic orbit starting from a uniform distribution of initial conditions. There are two straight nullclines in the phase portrait at = 0(vertical) and _ = 0(horizontal) lines. The intersection of these two nullclines is the system xed point at ; _ = (0;0). Figure 2.5(d) indicates the 3-dimensional view of the system phase portrait with respect to forcing displacement u = cos!t. In a periodic response we may use this 3D graph to count the number of forcing cycles in one period of the phase portrait cycle, which is the same as the number of fast time cycles in one slow time response period in the time trace. We can then obtain the frequency ratio of 28 : 1 from Figure 2.5(d). Notice that if the periodic orbit in the phase portrait is symmetric, then its ratio will have even numbers. For asymmetric phase portrait orbits, the ratio will have odd number. De nition 2 Let _ x =f(x) be a system of rst order ordinary di¤erential equation. Thex j nullcline is the set of points which satisfy f j (x 1 ;x 2 ;:::;x n ) = 0 and the intersection points of nullclines are equilibrium points of the system. From inspection of Figure 2.6(a) we may see the system time trace repeats itself in a period of around 0:3s; which refers to slow time of the system response. The fast time response of time trace is proportional to the periodic vertical forcing at normalized forcing parameters = 0:1710 and ! = 14:2784. From Figure 2.6(b) we can see that the system 34 Figure 2.6: Numerical integration results of single degree of freedom inverted pendulum under = 0:1710 and ! = 14:2784 with initial condition 0 = 0:1 and _ 0 = 0. (a)System time response. (b)Angular velocity time trace. (c)System phase portrait. (d)3D view of phase portrait with respect to periodic forcing. 35 angular velocity response also has the same repeated period around 0:3s: In the system phase portrait graph in Figure 2.6(c), we see a symmetric periodic orbit in a vector eld starting from uniform distribution of initial conditions. Figure 2.6(d) indicates 10 : 1 ratio from its 3D phase portrait. Figure 2.7: Numerical integration results of single degree of freedom inverted pendulum under normalized forcing amplitude = 0:2193 and normalized forcing frequency ! = 14:2784 with initial condition 0 = 0:1 and _ 0 = 0. (a)System time response. (b)Angular velocity time trace. (c)System phase portrait. (d)3D view of phase portrait with respect to periodic forcing. 36 From inspection of Figure 2.7(a) we may see the system time trace repeats itself in a period of around 0:25s; which refers to slow time of the system response. The fast time response of time trace is proportional to the periodic vertical forcing at normalized forcing parameters = 0:2193 and ! = 14:2784. From Figure 2.7(b) we can see that the system angular velocity response also has the same repeated period of around 0:25s: In the system phase portrait graph in Figure 2.7(c), we see an asymmetric periodic orbit in a vector eld starting from uniform distribution of initial conditions. Figure 2.7(d) indicates 7 : 1 ratio from its 3D phase portrait. From inspection of Figure 2.8(a) we may see the system time trace repeat itself in a period around 0:14s; which refers to slow time of the system response. The fast time response of time trace is proportional to the periodic vertical forcing at normalized forcing parameters = 0:400 and ! = 12:8506. From Figure 2.8(b) we can see that the system angular velocity response also has the same repeated period of around 0:14s: In the system phase portrait graph in Figure 2.8(c), we see a symmetric periodic orbit in a vector eld starting from uniform distribution of initial conditions. Figure 2.8(d) indicates 4 : 1 ratio from its 3D phase portrait. Since the system response is in general not periodic, more often we may get qua- siperiodic orbits for bounded solutions as shown in Figure 2.9. Figure 2.9 has forcing parameters = 0:200 and ! = 14:2784 which is very close to Figure 2.8. However, their system responses are total di¤erent. From time trace of Figure 2.9 we can only conclude that the system has bounded solution. Although system show symmetric pattern but we cannot obtain its repeatability as easy as in Figure 2.8. 37 Figure 2.8: Numerical integration results of single degree of freedom inverted pendulum undernormalizedforcingamplitude = 0:400andnormalizedforcingfrequency ! = 12:8506 with initial condition 0 = 0:56 and _ 0 = 0. (a)System time response. (b)Angular velocity time trace. (c)System phase portrait. (d)3D view of phase portrait with respect to periodic forcing. 38 Figure 2.9: Numerical integration results of single degree of freedom inverted pendulum undernormalizedforcingamplitude = 0:200andnormalizedforcingfrequency ! = 14:2784 with initial condition 0 = 0:1 and _ 0 = 0. (a)System time response. (b)Angular velocity time trace. (c)System phase portrait. (d)3D view of phase portrait with respect to periodic forcing. 39 Increasing the Forcing Amplitude Figure 2.10, 2.11 and 2.12 present the numerical results for the inverted pendulum undervaryingperiodicforcingamplitudes. Inthese3graphs, wevarytheforcingamplitude ; holding the normalized forcing frequency ! xed at ! = 14:2784. This normalized forcing frequency corresponds to the physical forcing frequency around 25Hz. All system initial conditions start from 1 0 = 0:1 and _ 1 0 = 0. Four system response representations are shown in a systematic way from left to right as: system response time trace, system phase portrait, 3D phase portrait, and power spectral density graph of system response, respectively. The normalized forcing amplitude increases from top to bottom and from Figure 2.10 to Figure 2.12. Recall Figure 2.4 in the previous chapter where this simulation data sets were presented in the sDOF inverted pendulum stability diagram shown in ! plane. This set of data points are referred to as numerical simulation 1 (Chen)indicated bythereddiamonddotsinFigure2.4. Asincreases, thelocationofthedatapointmoving toward the right side of Figure 2.4 horizontally with xed ! value. All data points were within the bounded region as predicted by the transition curves obtained from solving the linear Mathieus equation. From inspecting Figure 2.10 to Figure 2.12, we notice that as increases, the frequency ratio N r decreases. AsshownonthetopgraphinFigure2.10, under = 0:1019thetimetracehastwo distinct frequencies. One is referred to as the fast time response which is directly related to the forcing frequency; the other is referred to as slow time response which frequency is 54 times less than the forcing frequency. Its phase portrait has right moving vector eld when _ > 0; and left moving vector eld when _ < 0. Two nullclines are shown at the 40 horizontal line _ = 0 and the vertical line = 0 in the phase plane which indicate that the system has an equilibrium point located at _ ; = (0;0). The phase portrait shows one closed orbit, even with large numerical integration time. Its 3D phase portrait measured 54 : 1 frequency ratio where 54 forcing cycles can be counted from one phase portrait orbit. The power spectral density graph shows three distinct harmonic components, the largest peak in the low frequency region where the resolution of the graph is not ne enough to identify the exact location, the second peak is located at around 25Hz, and the third peak is located at around 50Hz. The low frequency peak is referred to as the slow time response of the system which may be measured from the time trace graph around 0:5Hz; the peak at around 25Hz is referred to as fast time system response which corresponds to the periodic forcing frequency. As we increase the value of , the systems slow time response becomes more rapid, resulting in a higher frequency in its largest peak in the power spectral density graph. If we inspect the evolution of system power spectral density graphs from increasing , we note that not only its largest peak moving toward the higher frequency region, its second and third peaks become wider in shape and eventually split into two new peaks as shown in = 0:1611. There were higher frequency peaks after the third peaks and they becamemoreobviouswhenincreased. Thesepeaksalsobecamewiderinshapeand nally split at higher values. When the frequency ratio N r has even value, the phase portrait exhibits a symmetric pattern, as shown in N r = 54;40;30;::;10;6;4; while the frequency ratio N r has odd value, the phase portrait exhibits an asymmetric pattern as shown in N r = 27;19;::;11;7;5. Figure 2.10 to Figure 2.12 does not show periodic solutions with every periodic. In fact, we may nd solutions in every period from N r 4 provided that 41 the numerical simulation tool has high enough resolution. Nonlinear Analysis of sDOF Model Using Perturbation Theory In order to understand the nonlinear e¤ect further in detail in the single degree of freedominvertedpendulummodeldescribedin(2:12),weappliedatwo-timingperturbation technique to the normalized model (2:56). Di¤erent than the direct numerical integration method that deal with time-varying model as described in previous sections, perturbation techniquecantransfertheoriginaltime-varyingequationofmotionintoasetofautonomous ordinary di¤erential equations. One of the advantages of applying perturbation technique is that large system response solutions can be obtained, system response wont be limited bythesmallangleassumptionasinthelinearcase. Also, thecostfornumericalsimulations on perturbation technique is usually less than direct numerical integration to the full time- varying model. We applied 2-timing perturbation technique by assigning t as slow system time and = !t as fast system time then making second time derivatives using chain rules: d 2 dt 2 = @ 2 @t 2 +2 ! @ 2 @t@ + ! 2 @ 2 @ 2 : (2.64) We then apply a uniform expansion using the normalized driving amplitude as a pertur- bation parameter in the following form: (t;;) = 0 (t;)+ 1 (t;)+ 2 2 (t;)+::::: (2.65) Substitute (2:64) and (2:65) into (2:56), 42 Figure 2.10: Numerical simulation results of periodic solutions from varying normalized forcing amplitude with xed normalized forcing frequency ! = 14:2784. Graphs from top to bottom show values increase from = 0:1019 to = 0:1240. From left to right showing four geometric representations of system response: time trace, phase portrait, 3D phase portrait, and PSD, respectively. 43 Figure 2.11: (Cont. of 2.10) Numerical simulation results of periodic solutions from vary- ing normalized forcing amplitude with xed normalized forcing frequency ! = 14:2784. Graphs from top to bottom show values increase from = 0:1293 to = 0:1531. From left to right showing four geometric representations of system response: time trace, phase portrait, 3D phase portrait, and PSD, respectively. 44 Figure 2.12: (Cont. of 2.11) Numerical simulation results of periodic solutions from vary- ing normalized forcing amplitude with xed normalized forcing frequency ! = 14:2784. Graphs from top to bottom show values increase from = 0:1611 to = 0:3346. From left to right showing four geometric representations of system response: time trace, phase portrait, 3D phase portrait, and PSD, respectively. 45 @ 2 @t 2 +2 ! @ 2 @t@ + ! 2 @ 2 @ 2 (t;;) | {z } LHS = 1+ ! 2 cos! n ! 2 n sin(t;;) | {z } RHS : (2.66) By expanding LHS(Left-hand-side) and RHS(Right-hand-side) of (2:66) using Taylor series expansion and then compare them with the same order, yields: O( 1 2 ) : @ 2 @ 2 0 = 0; (2.67) O( 1 ) : 2 @ 2 @t@ 0 +k @ 2 @ 2 1 =! 2 n [kcos! n sin 0 ]; (2.68) O(1) : @ 2 @t 2 0 +2k @ 2 @t@ 1 +k 2 @ 2 @ 2 2 =! 2 n sin 0 +k 2 cos! n ( 1 cos 0 ) ;(2.69) O() : @ 2 @t 2 1 +2k @ 2 @t@ 2 +k 2 @ 2 @ 2 3 = ! 2 n 1 cos 0 +k 2 cos! n ( 2 cos 0 1 2 2 1 sin 0 ) : (2.70) Equations(2:67);(2:68);(2:69)and(2:70)allful llthegoverningequationin(2:66). Notice that in above equations we introduced a scaling factor k which is de ned as normalized driving amplitude multiply by normalized driving frequency ! as k = ! (2.71) This k scale factor is related to the linear stability limit curve " 2 ! 2 = 2gL as shown in Figure 2.4. In our system k has value around 1:2 to be on the stability limit curve and stability condition in our system requires that k > 1:2. Since we are looking for bounded solutions (t) for (2:56), which we may get approximate solutions from solving the equations for 0 ; 1 ; 2 ;::: and then combine them using (2:65). The bounded solutions of (t) require that those equations of 0 ; 1 ; 2 ;::: are 46 also bounded. Using this bounded property we formulated the equations carrying out the integrations of individual equations to ful ll the original equation of motion. In the next sectionsweshowderivationsfromsolvingtheequationsfordi¤erentorderstogetthe nal approximate solution for (2:56): Solutions from Solving O( 1 2 ) Starting from the lowest order of equation, (2:67) is a homogeneous 2 nd order ordinary di¤erential equation which we may get solutions through direct integrations. We rst integrated the equation of order O( 1 2 ) with respect to , (2:67) becomes @ @ 0 (t;) =A(t) (2.72) where new function A(t) is being introduced from previous integration step which is not a function of . Making another integration to (2:72) with respect to yields 0 (t;) =A(t) +B(t) (2.73) where new function B(t) is being introduced from previous integration step. Since we were looking for solutions where 0 (t;) remain bounded as t; !1 , we must assign A(t) 0 which leads to 0 (t;) =B(t): (2.74) This result shows that the systems slow time response is synchronized with B(t) and the response is not a¤ected by the fast system time . So far, we can only conclude that B(t) is bounded in time and if we are looking for periodic solution of (t) then B(t) also needs 47 to be periodic in time. From carrying out integrations through higher order equations we may ndmoreconstraintsforB(t)toful lltheboundedconditions. Weshowedderivations from solving higher order equations in the following sections. Solutions from Solving O( 1 ) From (2:68) we can see that O( 1 ) equation coupled 0 and 1 in a 2 nd order ordinary di¤erential equation. The rst term in (2:68) can be calculated by making time derivatives to 0 (t;) with respect to t; yields @ 2 @t@ 0 (t;) = 0: (2.75) Substitute (2:75) and (2:74) into (2:68) yields @ 2 @ 2 1 =! 2 n cos! n sinB(t): (2.76) Integrate (2:76) with respect to we get @ @ 1 (t;) =! n sin! n sinB(t)+C(t) (2.77) wherenewfunctionC(t)isbeingintroducedfrompreviousintegrationstep. Makinganother integration to (2:77) with respect to yields 1 (t;) =cos! n sinB(t)+C(t)+D(t) (2.78) where new function D(t) is being introduced from previous integration step. Since our goal is to look for solutions that remain bounded as t; !1 , we may conclude that C(t) 0 48 is required to ful ll bounded conditions. Therefore, (2:78) becomes 1 (t;) =cos! n sinB(t)+D(t): (2.79) Since ! n =!t , we have 1 (t) =cos!tsinB(t)+D(t) (2.80) Therefore, we can conclude that 1 is a function of two new bounded function B(t) and D(t). Insearchingforperiodicsolutionsto(t)requiresthat 0 and 1 bothperiodicintime whichleadstosearchingforperiodicfunctionsofB(t)andD(t). Nextsectionwecontinued to carry out integrations to higher order equations in order to nd proper functions for B(t) and D(t). Solutions from Solving O(1) From (2:69) we can see that O(1) equation coupled 0 ; 1 and 2 in a 2 nd order ordinary di¤erential equation. We rearranged (2:69) as k 2 @ 2 @ 2 2 = @ 2 @t 2 0 2k @ 2 @t@ 1 +! 2 n sin 0 +k 2 cos! n ( 1 cos 0 ) : (2.81) Bycalculatingoneoftheintermediatetermsin(2:81)fromtakingtimederivativesto 0 (t;) with respect to t , we have @ 2 @t 2 0 (t;) = B(t): (2.82) Anotherintermediatetermin(2:81)canbecalculatedbytakingtimederivativesto @ @ 1 (t;) with respect to t 49 @ 2 @t@ 1 (t;) = _ B(t)! n sin! n cosB(t): (2.83) Substitute 0 and 1 equations in (2:74) and (2:80) into the 4 th term of (2:81), we have k 2 cos! n ( 1 cos 0 ) = k 2 cos! n (cos! n sinB(t)+D(t))cosB(t) = k 2 cos 2 ! n sinB(t)cosB(t) +k 2 D(t)cos! n cosB(t): (2.84) We then applied the following trigonometric identities to simplify (2:84) cos 2 ! n = 1 2 1 2 cos(2! n ); sinB(t)cosB(t) = 1 2 sin(2B(t)): Therefore (2:84) becomes k 2 cos! n ( 1 cos 0 ) = k 2 4 sin(2B(t))+ k 2 4 cos(2! n )sin(2B(t))+k 2 D(t)cos! n cosB(t): (2.85) Substitute (2:82),(2:83) and (2:85) into (2:81), we have 50 k 2 @ 2 @ 2 2 (t;) = B(t)2k _ B(t)! n sin! n cosB(t)+! 2 n sinB(t) k 2 ! 2 n 4 sin(2B(t))+ k 2 ! 2 n 4 cos(2! n )sin(2B(t)) +k 2 ! 2 n D(t)cos! n cosB(t) = B(t)! 2 n sinB(t)+ k 2 ! 2 n 4 sin(2B(t)) | {z } not function of 2k _ B(t)! n sin! n cosB(t)+ k 2 ! 2 n 4 cos(2! n )sin(2B(t)) +k 2 ! 2 n D(t)cos! n cosB(t): | {z } function of t; (2.86) From (2:86), bounded conditions lead to the following equation: B(t)! 2 n sinB(t)+ k 2 ! 4 n 4 sin(2B(t)) = 0; (2.87) D(t) = 0: (2.88) Equation (2:87) is an autonomous 2 nd order ordinary di¤erential equation governing the response ofB(t). Therefore, we may approximate the inverted pendulums response(t) by inserting the solution of (2:87) into the rst two terms of (2:65) (t) = 0 (t)+ 1 (t); = B(t)cos!tsinB(t): (2.89) We present numerical simulations of (2:87);(2:89) in later sections. 51 Solve B! 2 n sinB + k 2 ! 4 n 4 sin(2B) = 0 by Energy Equation In analyzing (2:87), we introduce the energy equation which we may obtain from multiplying both sides of (2:87) by _ B(t) as _ B(t) B(t)! 2 n _ B(t)sinB(t)+ k 2 ! 4 n 4 _ B(t)sin(2B(t)) = 0: (2.90) From (2:90) we get the following equation d dt 2 6 6 6 6 4 1 2 _ B(t) 2 | {z } Kinetic Energy T( _ B) + ! 2 n cosB(t) k 2 ! 4 n 8 cos(2B(t)) | {z } Potential Energy V(B) 3 7 7 7 7 5 = 0: (2.91) Equation (2:91) ful lls the energy conservation law as T +V =E =constant, (2.92) where system kinetic energy T takes the form T( _ B) = 1 2 _ B(t) 2 (2.93) and the system potential energy V takes the form V(B) = 2! 2 n cosB(t) k 2 ! 4 n 4 cos(2B(t)): (2.94) From (2:92) and (2:93) we can also get the following relation for _ B _ B(t) = p EV (2.95) where EV = T > 0: We present results by combining the potential energy graph to the solutions of (t) from solving B(t) numerically as shown in Figure 2.13. 52 Simulation Results Based on Perturbation Analysis Through the derivations from perturbation analysis we may transfer the origi- nal time-varying equation of motion as in (2:12) to describe the inverted pendulum under periodic forcing into a set of autonomous ordinary di¤erential equations as in (2:87) and (2:89). From solving (2:87) using numerical integration, and plugging the results of B(t) into (2:89), we may get the inverted pendulums response due to di¤erent forcing parame- ters. Figure 2.13 documents the inverted pendulums response with xed periodic forcing amplitude = 0:1233. Notice that we are now dealing with solving the autonomous system of B by treating the forcing parameters and ! as constants. The advantage of this system is that we can deal with large initial conditions and large system response. In Figure 2.13, each individual sub-graph indicates the condition of applying di¤erent forcing frequency and di¤erent initial conditions. Each sub-graph has four parts from top to bottom: system time response, system angular velocity response, system phase portrait with initial starting angle, and the potential energy graph at the starting initial conditions. The upper left graph in Figure 2.13 shows under normalized forcing frequency ! = 3:4 the system wont be stable in any initial conditions. The unstable behavior can be seen from its potential energy graph where the maximum potential occurs in its upright vertical position. There are three sub-graphs showing the conditions of ! = 10:2, the only stable case is when the system starts from initial conditions at 0 = 0:5236 and _ 0 = 0. These three system have the same system potential energy graph which indicates that the upright vertical position is the local minimum of the potential. As we increase the normalized forcing frequency !, we notice that the potential energy graphs change. In some cases the pendulum can be sta- 53 bilized near its inverted state even when the system doesnt have local minimum potential at its upright vertical position as in the case of ! = 20:4; 0 = 0:5236 and _ 0 = 0. Under proper conditions, the system will encounter bounded solutions as in ! = 17:4; 0 = 0:2356 and _ 0 = 0; insomecasesthesystemresponsesarenearperiodicasin ! = 20:4; 0 = 0:5236 and _ 0 = 0. Our simulation result showed that large periodic system response is possible even with large initial conditions. In the case of ! = 20:4; 0 = 0:8727 and _ 0 = 0, system started from around 50 initial angles and the pendulum was swinging in its inverted states with system response as large as 50 . Althoughourintentioninourresearchwasnttogetquantitativematchesbetween numerical simulations and experiments, our results showed that it is possible to get a close match from nd-tuned system parameters. Figure 2.14 shows one of the results of a comparison where the simulation followed closely to experiment. The top portion of the graphshowsthesystemtimetracewiththeinvertedpendulumsangleresponsewithrespect to time; bottom portion of the graph shows system phase portrait with the pendulums angle in the horizontal axis and the pendulums angular velocity in the vertical axis. Figure 2.14(a)shows the experimentaldatapoints of asingle degree of freedom inverted pendulum under vertical periodic forcing. The inverted pendulum was swinging in a large angle response with the largest swinging angle around 50 . Both time trace and phase portrait graphs indicate near periodic system behavior and the response showed very little damping e¤ect. Figure 2.14(b) shows the simulation result of applying perturbation technique to nonlinear single degree of freedom inverted pendulum model as in (2:56). The system has normalized forcing amplitude = 0:136 and normalized forcing frequency ! = 19:3 with 54 Figure 2.13: Numerical simulation based on perturbation analysis solutions. All simulation is under xed normalized driving amplitude = 0:1233. 55 initial conditions 0 = 45 and _ 0 = 0. 56 Figure 2.14: Quantitative comparison of the perturbation analysis simulation result to the experimental result. (a) ltered experimental time trace and phase portrait (b)simulation time trace and phase portrait. 57 Chapter 3 Experiments on the sDOF Inverted Pendulum In order to understand the inverted pendulum dynamics further, we conducted several experiments to compare in detail the results from measuring the responses of the real-life apparatus to the predictions from numerical simulations. In the following sections, we show the material and methods involved in our experiments, later, we present our ex- perimental results and then compare that to our numerical simulation results. 3.1 Material and Methods Driving Inverted Pendulum by Loud Speaker Following our theoretical model as described in (2:12), we conductedsDOF exper- iments with experimental apparatus as close to the assumptions in the theoretical model as possible. Recall our sDOF model as in Figure 2.1, our sDOF inverted pendulum model 58 consists of a point mass mounted on a weightless rigid rod which encounters periodic ver- tical forcing at its base. Since the ideal pendulum setup is not feasible experimentally, an alternate setup with light-weight rigid rod mounted with much heavier tip mass was used instead, as shown in Figure 3:1. The pivot point of the inverted pendulum was made by a commercially available high precision ball bearing. There are two major functions to this ball-bearing setup at its base: one is to minimize the rotational friction while the pendulum is in swinging motion. The other is to keep the pendulum moving on XY planer motion without any other 3 dimensional e¤ect. As in Figure 2.1, the pendulum consists of CNC machined rigid rod composed of Aluminum alloy with a heavier stainless steel made ball- bearing mounted on its tip. The full pendulum length is measured as L = 80:00:1mm, we measured its combined center of mass located at R = 53:70:1mm measured from the base with total weight of 3:80:1g. The base of the inverted pendulum was mounted at the ball-bearing pivot which was rmly attached to a stand mounted on top of the center cone of the speaker. The speaker was a 15in diameter commercially available subwoofer loud speaker. This loud speaker is capable of taking up to 140W RMS signals. The speaker core is composed of 2in thick voice coil with 40oz of magnet and its impedance rating is 4 . We applied a function generator to generate a sinusoidal wave and ampli ed it using a 100W PA power supply to drive the speaker. With this setup, it is capable of generating up to 50G of acceleration to drive the total weight of roughly 10g inverted pendulum apparatus to oscillate in a sinusoidal motion with forcing amplitude up to 10mm and forcing frequency up to 30Hz. The forcing amplitude and forcing frequency can be controlled by the function generator 59 Figure 3.1: Single degree of freedom inverted pendulum experimental setup. 60 and the ampli er. During each experiment, a pre-assigned forcing amplitude and forcing frequency was applied to drive the speaker and the pendulum base in a sinusoidal vibration while the pendulum was released. Under certain initial conditions, the pendulum stabilizes around its upright position swinging periodically with left and right motion without falling over. We then documented this pendulum motion with di¤erent forcing parameters and initial conditions in a systematic way and captured the pendulum motion using a high speed camera. Figure3.2: Preciseexperimentalmeasurementsfromimageacquisitionandimageprocessing using XS-3 high-speed camera, NI Vision Builder, MS Excel and Matlab. Inverted Pendulum Motion Using High Speed Camera In capturing the dynamics of an inverted pendulum, we are interested in mea- suring the input of the periodic forcing to the pendulum base and the resulting e¤ect of 61 the pendulum responses. By applying the optical measuring technique through high speed camera,wemeasuredthependulumsresponseswithhighprecision. Meanwhile,weavoided the problem of changing the sensitive inverted pendulums dynamics from applying tradi- tional mechanical measuring techniques which usually require some attachments from the sensors to the apparatus which often result in inducing system damping e¤ects. One IDT X-Stream XS-3 high speed and high image resolution camera was used for our experimental measurements. This high speed camera has xed image memory which we may trade in between of higher image resolution or longer recording time. We set the image resolution to 1260288 pixels with high camera frame rate which was recorded on 1G byte camera memory. The camera is capable of recording 3:6 seconds of data when setting the frame rate to 500 frames per second; 1:8 seconds of data when setting frame rate to 1000 frames per second. Sensing by applying the optical properties of camera would eliminate additional friction forces which would often be induced by mechanical sensors. However, sensing by the camera requires a large amount of computational power in image processing which results in slower processing time. As shown in Figure 3.2, there were four major steps of image processes involved in our experimental measurements of the inverted pendulum: (a)Image acquisition by XVision using XS-3 high speed camera, (b)Image processing by NI Labview and NI Vision Builder, (c) Data analysis by Microsoft Excel, and (d)Data ltering by Matlab Butterworth lter. The experimental setup was constructed as shown in Figure 3.3. An Agilent 33220A signal generator was used to generate sinusoidal forcing signal then a 100-Watt 62 PA ampli er was used to amplify the signal in order to drive the speaker so as to gener- ate periodic vertical forcing to the pendulum base. Notice that the whole speaker setup was mounted to a solid platform in order to minimize the table vibration. The pendu- lum is capable of adding additional pendulum rods through the ball-bearing joint, therefore higherdegreeoffreedomexperimentscanbeconductedfromsharingthesameexperimental equipment. Figure 3.3: Three dimensional view of the SDOF experimental setup. The motion of the pendulum was recorded by using an IDT X-Stream XS-3 high speed camera. To avoid image distortion form acquisition process, the camera view angle had been set to 0 with respect to the pendulum base and perpendicular to the pendulums swinging XY plane. The distance between the camera front lens and the pendulum base was 2m. The camera was set to acquire 8-bit resolution images in camera speed of 500 63 images per second and 1000 images per second. In order to enhance the image quality for better image analysis results, the rod of the pendulum and the pivot point on the base were painted with a reective white paint, while the rest of the apparatus was painted with a non-reective black coating to reduce glare. A ruler with a reective white paint on its front face and black coating on the tick marks was mounted in a xed location vertically onto the same pendulum moving plane. The distance between two tick marks on the ruler was 10mm away. As shown in Figure 3.1, the ruler had two important functions: one is to calibrate the images captured from high speed camera that we used the tick marks on the ruler to convert the pendulum size from image pixels to physical engineering unit in mm; the other function is to use the ruler to measure the pendulums angle where parallel to the ruler indicates 0 pendulum angle. Also, there was a DC powered truck headlight used to illuminate the apparatus in order to get better image quality. System Calibration and Image Processing by NI Vision Builder Since there were huge number of images captured by high speed camera taken from each experimental trial, it required some image processing steps to convert the im- age information into proper physical measurements. We applied NI Vision Builder image processing software to get our experimental results. NI Vision Builder can detect color(or brightness) di¤erence from an digital image, the shape of an object de ned by di¤erent color(or brightness) captured from the camera can then be distinguished. From scanning through a sequence of images from the same shape detecting algorithm in this software, we may get the time evolution of the shape changes which we de ned as motion of an object. Figure 3.4 shows a snap shot of a single image of a single degree of freedom inverted pendu- 64 lum in motion. Notice that the camera being mounted in a 90 angle in Figure 3.4 in order to get wider camera view from di¤erent aspect ratio to capture larger pendulums motion. The high speed camera provide black and white images with gray scale color and gradient brightness to distinguish the shape of the inverted pendulum. We rst applied NI Vision Builder to make more contrast to the image to enhance the image quality. The result is shown as in Figure 3.4, the inverted pendulum shows in white color which clearly separated from the black color background. Then we de ned an algorithm from NI Vision Builder to detect the shape of the pendulum base which we can assign a moving coordinate to follow the bases motion. We can then get the measurements of the pendulum bases motion in digital pixel coordinate as a function of image frames, later we converted this information into the time trace of the pendulums motion in engineering unit using mm=s. One of the zoom-in view of the pendulum bases motion is shown in Figure 3.5. Notice that we made use of a ruler with tick marks of 10mm to calibrate the image which can be seen from the top portion of Figure 3.4. Similarly, we detected the pendulum angle from capturing one of its straight edges in comparison of the vertical straight line provided from the ruler in the xed location. Those green boxes in Figure 3.4 indicate the regions of interest of the image process. Therefore, all pendulums responses and the forcing parameters captured fromhighspeedcameracanbemeasuredpreciselythroughaboveimageprocesstechniques. Data Analysis and Data Filtering Figure 3.5 shows the data points of the zoom-in view of one of the pendulums base motion. The Y signal in Figure 3.5 measured the pendulums vertical motion which indicated that this pendulum was under sinusoidal forcing with amplitude of 8:50:2mm 65 Figure 3.4: NI Vision Bulider software applied to a single degree of freedom experimental setup. and frequency of 25:0 0:1Hz. The X signal in Figure 3.5 measured the pendulums horizontal motion which indicated the pendulums base had small horizontal motion X = 0:00:1mm. Notice that Figure 3.5 indicated that the high speed camera had acquisition rate at 500 frames per second which had 50 images or data points shown within 0:1 second of time period. There were high frequencynoise induced from data acquisition and image process- ing steps. Experimental raw data cannot be properly analyzed without data ltering. In this paper, we applied Butterworth digital lter using Matlab to convert noisy raw signal intocleanerwaveforms. Imagedataacquisitionhadsamplingfrequencyat 500Hz, therefore its Nyquist frequency was at 250Hz. We applied 10 th order lowpass digital Butterworth lter to all of our experimental data with cuto¤ frequency being set to 15% of Nyquist 66 Figure 3.5: Zoom-in view of experimental data of inverted pendulum base motion. X denotes vertical periodic forcing displacement which has coordinate u(t) = "cos!t; Y denotes horiontal motion. 67 frequency. 3.2 Single Degree of Freedom System Experimental Results Below are graphs of our experimental results on sDOF inverted pendulum as de- scribed in Figure 3.1. In Figure 3.6, the pendulum was under sinusoidal forcing with ampli- tude 10:00:2mm and frequency 25:00:1Hz. We normalized the forcing amplitude and forcing frequency by the pendulums length and its natural frequency in order to compare withournumericalsimulationresults. Thesystemhasnormalizedamplitude " = 0:140:01 and normalized frequency ! = 13:30:2. Under 500Hz of high speed camera acquisition rate, we measured the pendulums response in high resolution. In part (a) of Figure 3.6 shows the pendulums vertical forcing with comparison to the pendulums response. It in- dicates that the input forcing was in periodic pattern near sinusoidal motion which resulted in a near periodic pendulum response. In part (b),(c),(d),(e),(f) of Figure 3.6 there are two parts of system response, the top graph indicates the raw data points while the bottom in- dicates the ltered data points by applying Butterworth lter. Notice from the graphs, the ltered data shown much cleaner results compare with raw data, especially in the system phase portrait graph in (e), the raw phase portrait cannot conclude much while the ltered phase portrait shows clear near periodic pattern. There were some damping e¤ect being detected from our measurements as in (c) of Figure 3.6. The system time trace in (c) shows that the inverted pendulums response had smaller and smaller amplitude in later timing. This damping e¤ect also implied that the pendulums response may eventually stabilized in its upright vertical position and the system was asymptotically stable with the damping 68 e¤ect. This result is di¤erent compare with our numerical simulation result in which we neglect damping through out. In (e), the system phase portrait shows symmetric pattern which also indicates the frequency ratio N r should be an even number. This pendulum was swinging in a special near periodic pattern which we characterized as 14 : 1 frequency ratio pattern. Thisexperimentalresultsmatchedwithournumericalsimulationresultsdescribed in earlier sections. In (f) of Figure 3.6 shows the power spectral density graph to extract the harmonic components of the pendulum response. In (f), it indicates that there was a slow time response greater than 1Hz and the next obvious peak happened at 25Hz which is corresponding to the sinusoidal forcing input and the fast system time. In Figure 3.7, the pendulum was under sinusoidal forcing with amplitude 8:0 0:3mm and frequency 25:00:1Hz. Similar to previous discussion, we normalized the forc- ing amplitude and forcing frequency by the pendulums length and its natural frequency in order to compare with our numerical simulation results. The system has normalized ampli- tude " = 0:110:01 and normalized frequency ! = 13:30:2. Under 500Hz of high speed camera acquisition rate, we measured the pendulums response in high resolution. In part (a) of Figure 3.7 shows the pendulums vertical forcing with comparison to the pendulums response. It indicates that the input forcing was in periodic pattern near sinusoidal motion which resulted in a near periodic pendulum response. In part (b),(c),(d),(e),(f) of Figure 3.7therearetwopartsofsystemresponse,thetopgraphindicatestherawdatapointswhile the bottom indicates the ltered data points by applying Butterworth lter. Notice from the graphs, the ltered data shown much cleaner results compare with raw data, especially in the system phase portrait graph in (e), the raw phase portrait cannot conclude much 69 Figure 3.6: Experiment on single degree of freedom pendulum under 10:00:2mm driving amplitude and 25Hz driving frequency. System has normalized amplitude " = 0:140:01 and normalized frequency ! = 13:30:2. 70 while the ltered phase portrait shows clear near periodic pattern. There were some damp- ing e¤ect being detected from our measurements as in (c) of Figure 3.7. The system time tracein(c)showsthattheinvertedpendulumsresponsehadsmallerandsmalleramplitude in later timing. This damping e¤ect also implied that the pendulums response may even- tually stabilized in its upright vertical position and the system was asymptotically stable with the damping e¤ect. This result is di¤erent compare with our numerical simulation result in which we neglect damping through out. In (e), the system phase portrait shows symmetric pattern which also indicates the frequency ratio N r should be an even number. This pendulum was swinging in a special near periodic pattern which we characterized as 30 : 1 frequency ratio pattern. This experimental results also matched with our numerical simulationresultsdescribedinearliersections. In(f)ofFigure3.7showsthepowerspectral densitygraphtoextracttheharmoniccomponentsofthependulumresponse. In(f),itindi- cates that there was a slow time response around 1Hz and the next obvious peak happened at 25Hz which is corresponding to the sinusoidal forcing input and the fast system time. In Figure 3.8, the pendulum was under sinusoidal forcing with amplitude 6:5 0:3mm and frequency 30:00:1Hz. Similar to previous discussion, we normalized the forc- ing amplitude and forcing frequency by the pendulums length and its natural frequency in order to compare with our numerical simulation results. The system has normalized ampli- tude " = 0:090:01 and normalized frequency ! = 15:90:2. Under 500Hz of high speed camera acquisition rate, we measured the pendulums response in high resolution. In part (a) of Figure 3.8 shows the pendulums vertical forcing with comparison to the pendulums response. It indicates that the input forcing was in periodic pattern near sinusoidal motion 71 Figure 3.7: Experiment on single degree of freedom pendulum under 8:00:3mm driving amplitude and 25Hz driving frequency. System has normalized amplitude " = 0:110:01 and normalized frequency ! = 13:30:2. 72 which resulted in a near periodic pendulum response. In part (b),(c),(d),(e),(f) of Figure 3.8therearetwopartsofsystemresponse,thetopgraphindicatestherawdatapointswhile the bottom indicates the ltered data points by applying Butterworth lter. Notice from the graphs, the ltered data shown much cleaner results compare with raw data, especially in the system phase portrait graph in (e), the raw phase portrait cannot conclude much while the ltered phase portrait shows clear near periodic pattern. There were some damp- ing e¤ect being detected from our measurements as in (c) of Figure 3.8. The system time tracein(c)showsthattheinvertedpendulumsresponsehadsmallerandsmalleramplitude in later timing. This damping e¤ect also implied that the pendulums response may even- tually stabilized in its upright vertical position and the system was asymptotically stable with the damping e¤ect. This result is di¤erent compare with our numerical simulation result in which we neglect damping through out. In (e), the system phase portrait shows symmetric pattern which also indicates the frequency ratio N r should be an even number. This pendulum was swinging in a special near periodic pattern which we characterized as 36 : 1 frequency ratio pattern. This experimental results also matched with our numerical simulationresultsdescribedinearliersections. In(f)ofFigure3.8showsthepowerspectral densitygraphtoextracttheharmoniccomponentsofthependulumresponse. In(f),itindi- cates that there was a slow time response around 1Hz and the next obvious peak happened at 30Hz which is corresponding to the sinusoidal forcing input and the fast system time. 73 Figure 3.8: Experiment on single degree of freedom pendulum under 6:50:3mm driving amplitude and 30Hz driving frequency. System has normalized amplitude " = 0:090:01 and normalized frequency ! = 15:90:2. 74 Chapter 4 Multiple Degree of Freedom Inverted Pendulum Model The multiple degree of freedom model supports di¤erent e¤ects than the single degree of freedom case. In this paper, we rst considered the general n degree of freedom inverted pendulum model under periodic vertical forcing at its base. In later sections, we validated sDOF and 2DOF models from assigning n = 1 and n = 2 to our n degree of freedom inverted pendulum model. 4.1 Formulate nDOF System Equation of Motion From La- grangian Approach Consider a nDOF inverted pendulum as shown in Figure 4.1 under planer motion with respect to a xed two dimensional coordinate XY. The inverted pendulum consists of n point masses which connect to each other by n rigid weightless rods. Assume each 75 joint is free to rotate on the XY plane without any friction or damping force. The base of the inverted pendulum is subject to vertical oscillatory forcing with vertical periodic displacement h(t) = "cos!t . For convenience we setup a sub-coordinate xy which has its origin (x 0 ;y 0 ) located on the inverted pendulums base. The n masses labelled from the rst joint counting from the base as m 1 ;m 2 ;:::;m k ;:::;m n with length of their rods label as l 1 ;l 2 ;:::;l k ;:::;l n . Each point mass has their coordinate with respect to the inertial frame written as (x 1 ;y 1 );(x 2 ;y 2 );:::;(x k ;y k );:::;(x n ;y n ) . The inverted pendulums motion can be described by the time evolution of the angles of the pendulum rods with respect to upright vertical position of the inertial frame and the pendulums angles are labelled as 1 (t); 2 (t); 3 (t);:::; n (t) respectively. The whole system is under the Earths gravitational acceleration g which is point downward. We can write the coordinate of each mass as x 1 = x 0 +l 1 sin 1 ; y 1 = y 0 +h(t)+l 1 cos 1 ; x 2 = x 0 +l 1 sin 1 +l 2 sin 2 ; y 2 = y 0 +h(t)+l 1 cos 1 +l 2 cos 2 ; ::: x n = x 0 +l 1 sin 1 +:::+l n sin n = n X i=1 l i sin i ; y n = y 0 +h(t)+l 1 cos 1 +:::+l n cos n =h(t)+ n X i=1 l i cos i : Taking time derivatives to the above equations we may get velocity components of each joint as 76 Figure 4.1: Geometric relationship of gereral case nDOF inverted pendulum. under vertical oscillatory forcing at its base 77 _ x 1 = l 1 _ 1 cos 1 ; _ y 1 = _ h(t)l 1 _ 1 sin 1 ; _ x 2 = l 1 _ 1 cos 1 +l 2 _ 2 cos 2 ; _ y 2 = _ h(t)l 1 _ 1 sin 1 l 2 _ 2 sin 2 ; ::: _ x n = l 1 _ 1 cos 1 +l 2 _ 2 cos 2 +:::+l n _ n cos n = n X i=1 l i _ i cos i ; _ y n = _ h(t)l 1 _ 1 sin 1 l 2 _ 2 sin 2 :::l n _ n sin n = _ h(t) n X i=1 l i _ i sin i : Therefore, the potential energy forn-DOF inverted pendulum can be expressed as V n = m 1 gy 1 +m 2 gy 2 +:::+m n gy n ; (4.1) = g 2 4 h(t) n X i=1 m i + n X i=1 n X j=1 m j l i cos i 3 5 : Kinetic energy for n-DOF inverted pendulum can be expressed as T n = m 1 2 _ x 2 1 + _ y 2 1 + m 2 2 _ x 2 2 + _ y 2 2 +:::+ m n 2 _ x 2 n + _ y 2 n ; (4.2) = 1 2 n X i=1 m i 2 6 6 6 4 _ h 2 (t)2 _ h(t) i P j=1 l i _ j sin j +2 i1 P j=1 l j l j+1 _ j _ j+1 cos( j j+1 )+2l i l 1 _ i _ 1 cos( i 1 )+ i P j=1 l 2 j _ 2 j 3 7 7 7 5 : The inverted pendulums equation of motion can then be derived by means of Lagrangesequationapplyingenergyconservationlaw. Forn-DOFsystemitrequiresnsets of Lagranges equations to describe the whole systems dynamic. We may write the k th Lagranges equation as 78 d dt @L n @ _ k @L n @ k = 0: (4.3) The system has n states 1 ; 2 ; 3 ;:::; n . In (4:3); k denotes the k th component of the states and _ k denotes the k th angular velocity component. The Lagrangian function de ned as L n = T n V n which is simply the system kinetic energy subtracts by the sys- tem potential energy. Lagranges equation describes the motion generated from the system kinetic energy will later become the system potential energy; the system potential energy will also transfer perfectly into kinetic energy then some motions being introduced. La- granges equation is bases on the energy conservation law where no energy being generated nor dissipated in a system. Unlike the Lagranges equation as in (2:2) described in earlier section that the vertical forcing being treated as nonconservative forces, our new system has the vertical forcing shown explicitly inside its kinetic energy. Dealing with conservative Lagrangesequationwillbene tfromsolvingonlyhomogeneousdi¤erentialequationswhich results in much lighter calculation cost. Refer to [62] and [61], Weibel and Baillieul intro- duced a more compact form of Lagrangian equation to their n degree of freedom normal pendulum( k angles are measured from downward vertical) under a rapidly forced cart on an inclined plane. We applied the compact form of Lagrangian equation to our n-DOF inverted pendulum as follow: L n (Q; _ Q;) = 1 2 _ Q T M (Q) _ Q+A(Q) T _ Q | {z } Kinetic Energy V (Q) |{z} : Potential Energy (4.4) The Lagrangian equation shown in (4:4) has no di¤erence compare with the La- grangian function calculated from our previous derivation on system kinetic energy and potential energy in (4:2) and (4:1) respectively. However, detail break down components of 79 system kinetic energy variation can be shown in more clear way in (4:4). As described in (4:4), the system kinetic energy is composed of two parts: the e¤ect of system inertial and the e¤ect of vertical periodic forcing. Here Q T = [ 1 2 3 ::: n ]2R n denotes the states of the system. _ Q is time derivative of Q denotes the angular velocity of the states. is the velocity component due to vertical periodic forcing and is described as (t) = _ h(t) ="!sin!t: (4.5) M (Q)2R nn is a state-dependent inertia tensor can be described in matrix form as: M (Q) = 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 M 11 () M 12 () : : M 1n () M 21 () : : : : : : : : M n1 () : : : M nn () 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 (4.6) where the components inside can be described as M ij () = 0 @ n X m k k=max(i;j) 1 A l i l j cos( i j ): (4.7) A(Q)2R n isknownasCorioliscouplingvectorwhichdescribestheverticalforcing a¤ect the kinetic energy of the system. A(Q) T = A 1 () A 2 () A 3 () ::: A n () (4.8) where the components inside can be described as A i () = n X k=1 m k ! l i sin i : (4.9) 80 The system potential energy V (Q)2R depends only on the states of system and can be described as V (Q) = n X i=1 " n X k=i m k ! gl i (cos i 1) # : (4.10) Although the general Lagrangian equation (4:4) is not necessarily conserved for arbitrary forcing input , the system can be considered conserved under periodic forc- ing. This had been pointed out by Weibel and Baillieul in [62],[61]. Our system only concern on the e¤ect of periodic forcing where (t +T) = (t) for some T > 0. In [62] and [61], Weibel and Baillieul introduced the concept of controller Hamiltonian H(Q;P;) corresponding to (4:4) via the Legendre transformation H(Q;P;) = P T _ Q L n where P T = P 1 P 2 P 3 ::: P n are canonical momenta de ned by P = @Ln @ _ Q . They con- cluded that for periodic forcing the controller HamiltonianH(Q;P;) can be averaged over one period of (t) to obtain the averaged Hamiltonian which results in a proper Hamil- tonian described as @H @t = 0 . The system can then be considered conserved thru averaging principle. Inlatersections, wewill rstverifythesingledegreeoffreedommodelbyapplying n = 1 to our n degree of freedom model and compare that to our earlier work as shown in (2:12). Two degree of freedom system will then being introduced after by applying n = 2 to our n degree of freedom model. Results on numerical simulation and experimentation of two degree of freedom system will be shown in later chapters and detailed comparison will then be performed. 81 4.2 Validate the Single Degree of Freedom Inverted Pendu- lum Model Apply n = 1 to (4:6);(4:8) and (4:10), we have Q T = [ 1 ]; (4.11) M (Q) = [M 11 ( 1 )] =m 1 l 2 1 ; (4.12) A(Q) = [A 1 ( 1 )] =m 1 l 1 sin 1 ; (4.13) V (Q) = m 1 gl 1 (cos 1 1): (4.14) Therefore, the single degree of freedom system Lagrangian equation can be ex- pressed as L 1 Q; _ Q; = 1 2 m 1 l 2 1 _ 2 1 ("!sin!t)(m 1 l 1 sin 1 ) _ 1 m 1 gl 1 (cos 1 1): (4.15) Derive @L 1 @ 1 ; @L 1 @ _ 1 and d dt @L 1 @ _ 1 as follow: @L 1 @ 1 = m 1 l 1 _ 1 ("!sin!t)cos 1 +m 1 gl 1 sin 1 ; (4.16) @L 1 @ _ 1 = m 1 l 2 1 _ 1 ("!sin!t)m 1 l 1 sin 1 ; (4.17) d dt @L 1 @ _ 1 = m 1 l 2 1 1 "! 2 cos!t m 1 l 1 sin 1 ("!sin!t)m 1 l 1 _ 1 cos 1 ; (4.18) The Lagranges equation becomes d dt @L 1 @ _ 1 @L 1 @ 1 =m 1 l 2 1 1 m 1 l 1 g+"! 2 cos!t sin 1 = 0: (4.19) We may divide both sides of (4:19) by m 1 l 2 1 to make its leading coe¢ cient equal to 1 yields 1 g+"! 2 cos!t 1 l 1 sin 1 = 0: (4.20) 82 The equation of motion in (4:20) is the same as in (2:12) derived from earlier sections. 4.3 Two Degree of Freedom Inverted Pendulum Model Apply n = 2 to (4:6);(4:8) and (4:10), we have Q T = 1 2 ; (4.21) M (Q) = 2 6 6 4 M 11 () M 12 () M 13 () M 14 () 3 7 7 5 = 2 6 6 4 (m 1 +m 2 )l 2 1 m 2 l 1 l 2 cos( 1 2 ) m 2 l 1 l 2 cos( 1 2 ) m 2 l 2 2 3 7 7 5 ; (4.22) A(Q) = A 1 A 2 = (m 1 +m 2 )l 1 sin 1 m 2 l 2 sin 2 ; (4.23) V (Q) = (m 1 +m 2 )gl 1 (cos 1 1)+m 2 gl 2 (cos 2 1): (4.24) Substitute (4:21);(4:22);(4:23) and (4:24) into (4:4), the Lagrangian function of 83 two degree of freedom inverted pendulum model can be expressed as L 2 Q; _ Q; = 1 2 _ 1 _ 2 2 6 6 4 M 11 () M 12 () M 13 () M 14 () 3 7 7 5 2 6 6 4 _ 1 _ 2 3 7 7 5 + A 1 A 2 2 6 6 4 _ 1 _ 2 3 7 7 5 V (Q); = 1 2 (m 1 +m 2 )l 2 1 _ 2 1 +m 2 l 1 l 2 _ 1 _ 2 cos( 1 2 )+ 1 2 m 2 l 2 2 _ 2 2 ("!sin!t) h (m 1 +m 2 )l 1 _ 1 sin 1 +m 2 l 2 _ 2 sin 2 i [(m 1 +m 2 )gl 1 (cos 1 1)+m 2 gl 2 (cos 2 1)]: (4.25) Derive @L 2 @ 1 ; @L 2 @ _ 1 and d dt @L 2 @ _ 1 as follow: @L 2 @ 1 = m 2 l 1 l 2 _ 1 _ 2 sin( 1 2 )("!sin!t)(m 1 +m 2 )l 1 _ 1 cos 1 +(m 1 +m 2 )gl 1 sin 1 ; (4.26) @L 2 @ _ 1 = (m 1 +m 2 )l 2 1 _ 1 +m 2 l 1 l 2 _ 2 cos( 1 2 ) ("!sin!t)(m 1 +m 2 )l 1 sin 1 ; (4.27) d dt @L 2 @ _ 1 = (m 1 +m 2 )l 2 1 1 +m 2 l 1 l 2 2 cos( 1 2 ) m 2 l 1 l 2 _ 2 _ 1 _ 2 sin( 1 2 ) "! 2 cos!t (m 1 +m 2 )l 1 sin 1 ("!sin!t)(m 1 +m 2 )l 1 _ 1 cos 1 : (4.28) The rst set of Lagranges equation becomes 84 d dt @L 2 @ _ 1 @L 2 @ 1 = (m 1 +m 2 )l 2 1 1 g+"! 2 cos!t (m 1 +m 2 )l 1 sin 1 +m 2 l 1 l 2 h 2 cos( 1 2 )+ _ 2 2 sin( 1 2 ) i = 0: (4.29) We may divide both sides of (4:29) by (m 1 +m 2 )l 2 1 to make its leading term 1 coe¢ cient equal to 1 yields 1 g+"! 2 cos!t 1 l 1 sin 1 + m 2 l 2 (m 1 +m 2 )l 1 h 2 cos( 1 2 )+ _ 2 2 sin( 1 2 ) i = 0: (4.30) Similarly, we may derive @L 2 @ 2 ; @L 2 @ _ 2 and d dt @L 2 @ _ 2 as follow: @L 2 @ 2 = m 2 l 1 l 2 _ 1 _ 2 sin( 1 2 )("!sin!t)m 2 l 2 _ 2 cos 2 +m 2 gl 2 sin 2 ; (4.31) @L 2 @ _ 2 = m 2 l 1 l 2 _ 1 cos( 1 2 )+m 2 l 2 2 _ 2 ("!sin!t)m 2 l 2 sin 2 ; (4.32) d dt @L 2 @ _ 2 = m 2 l 1 l 2 1 cos( 1 2 )m 2 l 1 l 2 _ 1 _ 1 _ 2 sin( 1 2 ) +m 2 l 2 2 2 "! 2 cos!t m 2 l 2 sin 2 ("!sin!t)m 2 l 2 _ 2 cos 2 : (4.33) The second set of Lagranges equation becomes d dt @L 2 @ _ 2 @L 2 @ 2 = m 2 l 1 l 2 1 cos( 1 2 )m 2 l 1 l 2 _ 2 1 sin( 1 2 ) +m 2 l 2 2 2 "! 2 cos!t m 2 l 2 sin 2 m 2 gl 2 sin 2 = 0: (4.34) 85 We may divide both sides of (4:34) by m 2 l 2 2 to make its leading term 2 coe¢ cient equal to 1 yields 2 g+"! 2 cos!t 1 l 2 sin 2 + l 1 l 2 h 1 cos( 1 2 ) _ 2 1 sin( 1 2 ) i = 0: (4.35) Therefore, two sets of equations (4:30) and (4:35) describe the dynamic of two degree of freedom inverted pendulum under vertical periodic forcing. However, these two setsofequationsarenotwritteninaconventionalformwhichhastwohighestorderstatesin thesameequation. Fortheconvenienceofapplyingconventionalnumericalsimulationtools, wemayfurthersimpli edtwosetsofequations (4:30)and (4:35)intogeneralformats. First setofsimpli edequationcanbeobtainedfrom (4:30)subtract (4:35)multipleby l 2 l 1 m 2 (m 1 +m 2 ) and then divide both sides of equation by its leading coe¢ cient 1 m 2 m 1 +m 2 cos 2 ( 1 2 ) , yields 1 g+"! 2 cos!t ( 11 sin 1 + 12 sin 2 )+ 13 _ 2 1 + 14 _ 2 2 = 0: (4.36) The coe¢ cients inside (4:36) are 11 = 1 l 1 h 1 m 2 m 1 +m 2 cos 2 ( 1 2 ) i; (4.37) 12 = 1 l 1 m 2 (m 1 +m 2 ) cos( 1 2 ) h 1 m 2 m 1 +m 2 cos 2 ( 1 2 ) i; (4.38) 13 = m 2 (m 1 +m 2 ) cos( 1 2 )sin( 1 2 ) h 1 m 2 m 1 +m 2 cos 2 ( 1 2 ) i ; (4.39) 14 = l 2 l 1 m 2 (m 1 +m 2 ) sin( 1 2 ) h 1 m 2 m 1 +m 2 cos 2 ( 1 2 ) i: (4.40) 86 Secondsetofsimpli edequationcanbeobtainedfrom(4:35)subtract(4:30)multi- pleby l 1 l 2 andthendividebothsidesofequationbyitsleadingcoe¢ cient1 m 2 m 1 +m 2 cos 2 ( 1 2 ) , yields 2 g+"! 2 cos!t ( 21 sin 1 + 22 sin 2 )+ 23 _ 2 1 + 24 _ 2 2 = 0: (4.41) The coe¢ cients inside (4:41) are 21 = 1 l 2 h 1 m 2 m 1 +m 2 cos 2 ( 1 2 ) i; (4.42) 22 = 1 l 2 cos( 1 2 ) h 1 m 2 m 1 +m 2 cos 2 ( 1 2 ) i; (4.43) 23 = l 1 l 2 sin( 1 2 ) h 1 m 2 m 1 +m 2 cos 2 ( 1 2 ) i; (4.44) 24 = m 2 (m 1 +m 2 ) cos( 1 2 )sin( 1 2 ) h 1 m 2 m 1 +m 2 cos 2 ( 1 2 ) i : (4.45) Equations(4:36)and(4:41)describethedynamicoftwodegreeoffreedominverted pendulum under vertical periodic forcing at its base. 4.4 Coe¢ cients of 2DOF Inverted Pendulum Equation of Motion Theequationofmotionof2DOFinvertedpendulumunderperiodicverticalforcing is described by (4:36) and (4:41). These two sets of 2 nd order nonlinear ordinary di¤erential equations has coupled time-varying coe¢ cients 11 ; 12 ; 21 ; 22 due to periodic forcing. In addition, there are kinetic energy terms being induced through the interaction of two 87 pendulumsmotiontoeachother,theircoe¢ cientsare 13 ; 14 ; 23 ;and 24 . Noticethatall coe¢ cientsarestate-dependentandtheyalldependonthedi¤erencebetweentwopendulum angles 1 2 . Due to sin( 1 2 ) term in the coe¢ cients of 13 ; 14 ; 23 ; and 24 , their values are always small since the angle 1 2 is always small when pendulum is stabilized in its inverted states. We may also notice that all coe¢ cients in the system equation of motion share the same denominator as 1 m 2 m 1 +m 2 cos 2 ( 1 2 ). This denominator has value smaller than unity at all conditions which is also state-dependent 0< m 1 m 1 +m 2 1 m 2 m 1 +m 2 cos 2 ( 1 2 ) 1: (4.46) Since all parameters has positive values m 1 ;m 2 ;l 1 ;l 2 > 0 , the coe¢ cient 11 in the rst set of equation of motion (4:36) is always larger than 12 because 1 l 1 > 1 l 1 m 2 (m 1 +m 2 ) cos( 1 2 ) (4.47) due to the facts that m 2 (m 1 +m 2 ) < 1 and1 < cos( 1 2 ) < 1. The coe¢ cient 11 always has positive value. From (4:37) and (4:38), we can also get the following relationship: m 2 (m 1 +m 2 ) 11 12 m 2 (m 1 +m 2 ) 11 (4.48) which shows the limit of 12 values. Similarly, the coe¢ cient 21 in the second set of equation of motion (4:41) is always larger than 22 . From (4:39) and (4:45) we can get the following relationship: 13 = 24 (4.49) 88 which states that 13 has the same quantity as 24 but has reverse sign. Notice that if we set m 2 ! 0, then 12 ; 13 ; 14 approach zero value which results in the same equation of motion as in single degree of freedom model (2:13): This is the special case of the 2DOF inverted pendulum. Two degree of freedom pendulum has two distinct modes of vibration when it stabilized in its inverted state under vertical periodic forcing . Figure 4.2 shows two modes of vibration. The graph on the left shows two pendulums swinging in the same direction where _ 1 and _ 2 has the same sign, we named this kind of vibration as mode 1 type. On the other hand, the graph on the right shows two pendulums swinging in the opposite directions where _ 1 and _ 2 have di¤erent signs, we named this kind of vibration as mode 2 type. Although these two types of vibration modes have di¤erent pendulum responses they share the same equation of motion sets as in (4:36) and (4:41). 4.5 NumericalIntegrationonNonlinear2DOFInvertedPen- dulum Model Introduce new system states as x 1 = 1 ;x 2 = 2 ;x 3 = _ 1 ;x 4 = _ 2 and two additional states u = cos!t and v = sin!t that adapt to periodic forcing terms, then nonlinear 2DOF system equation of motion (4:36) and (4:41) becomes 89 Figure 4.2: Simpli ed two degree of freedom inverted pendulum model with two point masses mounted on two rigid rods. Two modes of vibration shown on graph: Mode 1(left) and Mode 2(right). _ x 1 = x 3 ; _ x 2 = x 4 ; _ x 3 = g+"! 2 u ( 11 sinx 1 + 12 sinx 2 ) 13 x 2 3 + 14 x 2 4 ; _ x 4 = g+"! 2 u ( 21 sinx 1 + 22 sinx 2 ) 23 x 2 3 + 24 x 2 4 ; _ u = u(1u 2 v 2 )!v; _ v = v(1u 2 v 2 )+!u: (4.50) All coe¢ cients inside (4:50) are the same de nitions as in equations from (4:37) to (4:45). We may apply numerical integration technique through parameter continuation code using 4 th order Runge-Kutta scheme to (4:50). Two simulation results are shown in Figure 4.3 90 and Figure 4.4. Both simulation results were obtained from applying the same system parameters: m 1 = 0:0038;m 2 = 0:0019;l 1 = 0:0537;l 2 = 0:0251. Figure 4.3 shows that systemisundermode1vibrationwherethe rstpendulumisswingingatthesamedirection as the second pendulum. Figure 4.4 shows that system is under mode 2 vibration where the rst pendulum is swinging at opposite directions compare to the second pendulum. Both simulationresultscamefromthesameverticalperiodicforcingtothependulumbasewhich has periodic forcing amplitude " = 0:00593 and forcing frequency ! = 250rad=s. The only di¤erence in getting these two distinct responses from the same equation of motion sets are applying di¤erent initial conditions to (4:50). Mode 1 response or Figure 4.3 has initial conditions 1 0 = 0:05; _ 1 0 = 0 and 2 0 = 0:1; _ 2 0 = 0; Mode 2 response or Figure 4.4 has initial conditions 1 0 = 0:1; _ 1 0 = 0 and 2 0 =0:101; _ 2 0 = 0. As shown in Figure 4.3, the inverted pendulum has bounded response. Its time trace represents near periodic behavior and the rst pendulum swings at the same direction asthesecondpendulumwhichindicatedmode1vibration. Bothofitsphaseportraitsshow bounded orbits which we may conclude as quasiperiodic behavior. Unlike single degree of freedom inverted pendulum case whose phase portrait has zero nullclines in both vertical and horizontal directions, the 2DOF inverted pendulum has non-zero vertical nullcline if it starts from non-zero initial conditions. The starting vertical nullcline location in 1 _ 1 plane depends on the initial condition of the second pendulum while the starting vertical nullclinelocationin 2 _ 2 planedependsontheinitialconditionsofthe rstpendulum. In fact, thoseverticalnullclines obtainedfromthenonlinearmodelwere notperfectlystraight, their curvatures depended on the coe¢ cients 13 ; 14 ; 23 ; 24 . The 3D phase portrait in 91 Figure 4.3: 2DOF inverted pendulum under mode 1 vibration with vertical periodic forcing amplitude " = 0:00593 and forcing frequency ! = 250rad=s. Initial conditions are 1 0 = 0:05; _ 1 0 = 0 and 2 0 = 0:1; _ 2 0 = 0 . 92 Figure 4.3 shown that their orbits were close to 4 : 1 periodic ratio. The power spectral densitygraphshownthatthesystemresponsehadalreadypastperiod-doublingbifurcation, there were multiple harmonic components higher than the forcing frequency which was ! = 250rad=s 39:8Hz. Compare to Figure 4.3, when applying di¤erent initial conditions to the same 2DOF inverted pendulum model with the same forcing as in (4:50) we can get totally di¤erent pendulum response as in Figure 4.4 which indicates the 2DOF inverted pendulum is now under mode 2 vibration. In Figure 4.4 the pendulum response can only be concluded as bounded, the periodic e¤ect is not obvious. Its time trace shown the rst pendulum swings at di¤erent directions compare with the second pendulum which indicated mode 2 vibration. Both of its phase portraits and its 3D phase portrait shown bounded behavior. Notice that although two power spectral density graphs in Figure 4.3 and Figure 4.4 showed similar patterns. However, in Figure 4.3 it showed the dominated harmonic component located at around 10Hz; while in Figure 4.4 the dominated harmonic component had been shifted to lower frequency at around 2Hz. None of our simulation results shown clear periodic behavior in 2DOF inverted pendulum using full nonlinear model as described in (4:50). This can be explained as the induced kinetic terms a¤ected the steadiness of their phase portrait orbits and caused them to uctuate regardless that their coe¢ cients 13 ; 14 ; 23 ; 24 are always small. 93 Figure4.4: Whenapplyingdi¤erentinitialconditionstothesame2DOFinvertedpendulum model we can get mode 2 vibration with the same vertical periodic forcing amplitude " = 0:00593 and forcing frequency ! = 250rad=s. Initial conditions are 1 0 = 0:1; _ 1 0 = 0 and 2 0 =0:101; _ 2 0 = 0. 94 4.6 Small Angle Assumption to Obtain Linear 2DOF In- verted Pendulum Model To obtain clear periodic orbits in a 2DOF case, we may linearize its nonlinear model by applying small angle assumptions to pendulum angles. Under the condition of small 1 and 2 angles as 1 ! 0 and 2 ! 0, we may get the following approximations: sin 1 ! 1 ; sin 2 ! 2 ; cos 1 ! 1; cos 2 ! 1; (4.51) and sin( 1 2 ) ! ( 1 2 ); cos( 1 2 ) ! 1; cos 2 ( 1 2 ) ! 1: (4.52) Therefore, we may obtain the simpli ed system equation of motion as 1 g+"! 2 cos!t ( 11 1 + 12 2 )+ 13 _ 2 1 + 14 _ 2 2 = 0; (4.53) 2 g+"! 2 cos!t ( 21 1 + 22 2 )+ 23 _ 2 1 + 24 _ 2 2 = 0: (4.54) where the new coe¢ cients can be simpli ed as 95 11 ! m 1 +m 2 m 1 l 1 ; (4.55) 12 ! m 2 m 1 l 1 ; (4.56) 13 ! m 2 m 1 ( 1 2 ); (4.57) 14 ! m 2 l 2 m 1 l 1 ( 1 2 ); (4.58) and 21 ! m 1 +m 2 m 1 l 2 ; (4.59) 22 ! m 1 +m 2 m 1 l 2 ; (4.60) 23 ! (m 1 +m 2 )l 1 m 1 l 2 ( 1 2 ); (4.61) 24 ! m 2 m 1 ( 1 2 ): (4.62) Since 1 and 2 are both small, the coe¢ cients 13 ; 14 ; 23 ; 24 are also small. We may then neglect those induced kinetic terms due to 13 ; 14 ; 23 ; 24 ! 0 to obtain linear system equation of motion as: 1 g+"! 2 cos!t ( 11 1 + 12 2 ) = 0 (4.63) 2 g+"! 2 cos!t ( 21 1 + 22 2 ) = 0 (4.64) with new system coe¢ cients 96 11 = m 1 +m 2 m 1 l 1 ; (4.65) 12 = m 2 m 1 l 1 ; (4.66) 21 = m 1 +m 2 m 1 l 2 ; (4.67) 22 = m 1 +m 2 m 1 l 2 : (4.68) Apply change of variables ofx 1 = 1 , x 2 = 2 , x 2 = _ 1 , x 4 = _ 2 , the linear system equation of motion can be described in matrix form as 2 6 6 6 6 6 6 6 6 6 6 4 _ x 1 _ x 2 _ x 3 _ x 4 3 7 7 7 7 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 6 6 6 6 4 0 0 1 0 0 0 0 1 (t) 11 (t) 12 0 0 (t) 21 (t) 22 0 0 3 7 7 7 7 7 7 7 7 7 7 5 2 6 6 6 6 6 6 6 6 6 6 4 x 1 x 2 x 3 x 4 3 7 7 7 7 7 7 7 7 7 7 5 : (4.69) Equation (4:69) takes the form _ X = A(t)X which is a rst order ordinary di¤erential equation with time-periodic coe¢ cient. The time-periodic property is coming from the vertical periodic forcing h(t) = "! 2 cos!t which is included inside (t) as (t) = g + "! 2 cos!t . Equation (4:69) is similar to (2:50) which indicates that Floquet theory and Hills equation can be applied to check the stability of this 2DOF linear time-periodic system. It is possible to obtain clear periodic orbits from simulating the linear 2DOF invertedpendulummodelin(4:69)sinceoneofthee¤ectofcausinguctuationsinnonlinear model has been removed from making the coe¢ cients of _ 2 1 and _ 2 2 equal to zero. Since there are two modes of vibration in the 2DOF inverted pendulum response, wemayapplythesametechniqueasinsolvingthenonlinearmodel(4:50)numericallywhich requires to assign two di¤erent types of initial conditions: assign 1 0 and 2 0 the same sign 97 and _ 1 0 = _ 2 0 = 0 to get mode 1 vibration; assign 2 0 reverse sign to 1 0 and _ 1 0 = _ 2 0 = 0 to get mode 2 vibration. Notice that in mode 2 case when we reverse sign of two initial angles, it also requires that coe¢ cients 21 ; 22 in second set of equation of motion (4:64) reverse their signs as shown in (4:70) and (4:71). 21 = m 1 +m 2 m 1 l 2 ; (4.70) 22 = m 1 +m 2 m 1 l 2 : (4.71) Our goal in numerical simulation on 2DOF linear inverted pendulum model is to characterize the type of periodic solutions similarly to the single degree of freedom case presented in earlier chapter of this paper. We applied the same parameter continuation simulationcodeasbeingusedinsDOFcasebyintegratingtheordinarydi¤erentialequations numerically through 4 th order Runge-Kutta scheme. The simulation results of two modes of vibration on 2DOF linear inverted pendulum under periodic forcing is presented in the following sections. Similar to sDOF case, four di¤erent kinds of geometric realization of system response are being used to present our data: system time traces, system phase portrait, three dimensional view of system phase portrait compare to input forcing, and systems power spectral density diagram. SimulationResultsonLinear2DOFInvertedPendulumModelunderMode 1 Vibration Since our goal in numerical simulation is to characterize di¤erent types of periodic solutions due to di¤erent periodic forcing, we applied a set of generic coe¢ cients as 11 = 98 1:5; 12 = 0:6; 21 = 0:6; 22 = 0:6 to our 2DOF linear model as described in (4:63) and (4:64). All simulations were under the same initial conditions as 1 0 = 0:1; _ 1 0 = 0 and 2 0 = 0:3; _ 2 0 = 0. Two periodic solution results are shown in Figure 4.5 and Figure 4.6. The system responses indicate that both system were under mode 1 vibration. We used the same four geometric representations: time trace, phase portrait, 3D phase portrait, and power spectral density graph, to present our simulation results. The left hand side graph setindicatesthe rstdegreeoffreedompendulumsresponsewhilerighthandsidegraphset indicates the second degree of freedom pendulums response. In Figure 4.5, with periodic input forcing amplitude " = 0:2138 and forcing frequency ! = 200 it shows that pendulum was under 4 : 1 frequency ratio. Its phase portraits represent symmetric periodic orbits withtheirverticalnullclineslocatedinnon-zeropositions. Thepowerspectraldensitygraph extracts multiple peaks in various frequencies locations and the largest peak happened in the lowest frequency at around 10Hz. In Figure 4.6, with periodic input forcing amplitude " = 0:03025 and forcing frequency ! = 200 it shows that pendulum was under 35 : 1 frequency ratio. Its phase portraits have asymmetric periodic orbits and their vertical nullclines located in non-zero positions. The power spectral density graph extracts 3 clear peaks with the dominated harmonic component located below 1Hz denotes as system slow time and second peak corresponding to the periodic forcing around 32Hz denotes as system fast time. The third peak has frequency about twice as large as the second peak which located within the range of 60Hzs 70Hz: Compare with periodic response as shown in Figure 4.6, under slightly di¤erent 99 Figure 4.5: 2DOF inverted pendulum under mode 1 vibration. The periodic input forcing amplitude was " = 0:2138 and the forcing frequency was ! = 200. 100 Figure 4.6: 2DOF inverted pendulum under mode 1 vibration. The periodic input forcing amplitude was " = 0:03025 and the forcing frequency was ! = 200. 101 periodic input forcing amplitude" = 0:03023, the pendulums response doesnt repeat itself exactly as in previous periodic case. Instead, the system showed quasiperiodic behavior as in Figure 4.7. In part (a) of Figure 4.7, we presented time trace to a longer period of 10 seconds in order to see long term signal drift. Its phase portraits show that the system response still remain bounded but the orbits in the phase portrait occupy more area inside. Their 3D phase portraits also occupy more area in the graph. The power spectral density graph is almost identical to the previous case as shown in Figure 4.6 which extracts 3 clear peaks with the dominated low frequency harmonic component located below 1Hz denotes as system slow time and second peak corresponding to the periodic forcing around 32Hz denotes as system fast time. The third peak has frequency about twice as large as the second peak which located within the range of 60Hzs 70Hz: SimulationResultsonLinear2DOFInvertedPendulumModelunderMode 2 Vibration Similar to the simulation results on linear 2DOF inverted pendulum model under mode 1 vibration, we applied a set of generic coe¢ cients to simulate mode 2 vibration as 11 = 1:5; 12 = 0:6; 21 = 0:6; 22 = 0:6 to our 2DOF linear model as described in (4:63)and(4:64). Allsimulationswereunderthesameinitialconditionsas 1 0 = 0:1; _ 1 0 = 0 and 2 0 =0:05; _ 2 0 = 0. Four periodic solution results are shown in Figure 4.8 , Figure 4.9, Figure 4.9, and Figure 4.11. The system responses indicate that they were under mode 2 vibration. We used the same four geometric representations: time trace, phase portrait, 3Dphaseportrait, andpowerspectraldensitygraph, topresentoursimulationresults. The left hand side graph set indicates the rst degree of freedom pendulums response while 102 Figure 4.7: 2DOF inverted pendulum under mode 1 vibration. The periodic input forcing amplitude was " = 0:03023 and the forcing frequency was ! = 200. The inverted pendulum shows quasiperiodic response. 103 right hand side graph set indicates the second degree of freedom pendulums response. In Figure 4.8, with periodic input forcing amplitude " = 0:2138 and forcing frequency! = 200 it shows that pendulum was under 4 : 1 frequency ratio. This is the same result as shown in mode 1 case as in Figure 4.5. Its phase portraits represent symmetric periodic orbits with their vertical nullclines located in non-zero positions. Notice that the second degree of freedom phase portrait has its non-zero nullcline moved to the opposite location compare to mode 1 case as in Figure 4.5. The power spectral density graph extracts multiple peaks in various frequencies locations and the largest peak happened in the lowest frequency at around 10Hz. In Figure 4.9, with periodic input forcing amplitude " = 0:1514 and forcing fre- quency! = 200itshowsthatpendulumwasunder6 : 1frequencyratio. Itsphaseportraits have symmetric periodic orbits and their vertical nullclines located in non-zero positions. Thepowerspectraldensitygraphextractsmultiplepeakswiththedominatedlowfrequency harmonic component located around 6Hz denotes as system slow time. In Figure 4.10, with periodic input forcing amplitude " = 0:11637 and forcing frequency ! = 200 it shows that pendulum was under 8 : 1 frequency ratio. Its phase portraits have symmetric periodic orbits and their vertical nullclines located in non-zero positions. The power spectral density graph extracts multiple peaks similar to 6 : 1 case with the dominated low frequency harmonic component moved to even lower frequency around 4Hz denotes as system slow time. In Figure 4.11, with periodic input forcing amplitude " = 0:05037 and forcing frequency ! = 200 it shows that pendulum was under 18 : 1 frequency ratio. Its phase 104 Figure 4.8: 2DOF inverted pendulum under mode 2 vibration. The periodic input forcing amplitude was " = 0:2138 and the forcing frequency was ! = 200. 105 Figure 4.9: 2DOF inverted pendulum under mode 2 vibration. The periodic input forcing amplitude was " = 0:1514 and the forcing frequency was ! = 200. 106 Figure 4.10: 2DOF inverted pendulum under mode 2 vibration. The periodic input forcing amplitude was " = 0:11637 and the forcing frequency was ! = 200. 107 portraits have symmetric periodic orbits and their vertical nullclines located in non-zero positions. The power spectral density graph extracts 3 clear peaks with the dominated low frequency harmonic component moved to even lower frequency around 1:5Hz denotes as system slow time. The second peaks located at around 32Hz corresponding to periodic forcing frequency and the system fast time. The third peak has frequency about twice as large as the second peak which located within the range of 60Hzs 70Hz: Compare with periodic response as shown in Figure 4.11, under slightly di¤erent periodic input forcing amplitude " = 0:0500, the pendulums response doesnt repeat itself exactly as in previous periodic case. Instead, the system showed quasiperiodic behavior as in Figure 4.12. In part (a) of Figure 4.12, we presented time trace to a longer period of 10 seconds in order to see long term signal drift. Its phase portraits show that the system response still remain bounded but the orbits in the phase portrait occupy more area inside. Their 3D phase portraits also occupy more area in the graph. The power spectral density graph is almost identical to the previous case as shown in Figure 4.11 which extracts 3 clear peaks with the dominated low frequency harmonic component located around 1:5Hz denotes as system slow time and second peak corresponding to the periodic forcing around 32Hz denotes as system fast time. The third peak has frequency about twice as large as the second peak which located within the range of 60Hzs 70Hz: Compare with periodic response as shown in Figure 4.8, under di¤erent periodic input forcing amplitude " = 0:2600, the pendulums response still performed near periodic behavior as in Figure 4.13. In order to see long term signal drift we presented time trace to a longer period of 20 seconds as shown in part (a) of Figure 4.13. Its phase portraits 108 Figure 4.11: 2DOF inverted pendulum under mode 2 vibration. The periodic input forcing amplitude was " = 0:05037 and the forcing frequency was ! = 200. 109 Figure 4.12: 2DOF inverted pendulum under mode 2 vibration. The periodic input forcing amplitude was " = 0:05 and the forcing frequency was ! = 200. The inverted pendulum shows quasiperiodic response. 110 showedwiderperiodicorbitsin 4 : 1frequencyratioasinpart(c),(d)ofFigure4.13andthe phase portrait orbits were symmetric with non-zero vertical nullclines. The power spectral densitygraphsarealmostidenticaltothepreviouscaseasshowninFigure4.8whichextract multiple peaks in various frequencies locations and the largest peak happened in the lowest frequency at around 10Hz. Numerical simulation result comparison From analyzing the response of linear 2DOF inverted pendulum model we found two modes of vibration from di¤erent initial conditions. However, in searching for periodic solutions inside the stable region around pendulums inverted state, our simulation results showed similar patterns from comparing these two modes. Therefore we choose only one mode of simulation results in presenting graph comparison from di¤erent forcing amplitude while keeping forcing frequency xed. Figure 4.14 presents the periodic solutions of system phase portrait comparison of the same 2DOF linear model under increasing forcing am- plitudes. Notice that in Figure 4.14 we applied the same periodic forcing frequency ! = 200rad=s which is equivalent to about 32Hz throughout. From upper left to lower right of Figure 4.14, we applied ascending input periodic forcing amplitudes from " = 0:04015 to " = 0:259. Their phase portraits show di¤erent periodic orbits and their frequency ratios have descending orders from increasing ". Notice that we dont list all of the periodic solu- tions, only 9 graphs being chosen that we can display them in a systematic way. The lowest frequency ratio that we got from simulation results is 4 : 1. CorrespondingtothephaseportraitcomparisongraphsasinFigure4.14,inFigure 4.15 we show their power spectral density graphs in relative graph locations according to 111 Figure 4.13: 2DOF inverted pendulum under mode 2 vibration. The periodic input forcing amplitude was " = 0:26 and the forcing frequency was ! = 200. The inverted pendulum shows near periodic response. 112 Figure 4.14: Phase portrait comparison of periodic solutions on numerical simulation using linear2DOFinvertedpendulummodel. Fromupperlefttolowerright,theforcingamplitude has ascending values which result in descending frequency ratio. 113 the same ascending order of periodic forcing amplitudes from " = 0:04015 to " = 0:259. From upper left graph of Figure 4.15, the system is under low periodic forcing amplitude at " = 0:04015 and its power spectral density graph shows 3 clear peaks. The largest peak happened at the lowest frequency location around 3Hz which is referred as slow time of the system. The second peak has frequency around 30Hz which follows the periodic forcing frequency which is also referred as fast time of system. The third peak has about double frequency compare to the second one. As we increased the periodic forcing amplitude ", the slow time response moved to slightly higher frequencies and the second and third peaks eventually split into two new frequencies in their neighborhood. Also notice that there were higher frequency peaks (the fourth peak) grew in higher periodic forcing amplitude case as in " = 0:09476. This fourth peak has frequency about 3 times larger than the second peak located at around 100Hz. In Figure 4.15, we can also see that the fourth peak split into two frequencies in the higher forcing amplitude cases. 114 Figure 4.15: Power spectral density comparison of periodic solutions on numerical simu- lation using linear 2DOF inverted pendulum model. From upper left to lower right, the forcing amplitude has ascending values which result in descending frequency ratio. 115 Chapter 5 Experiments on the 2DOF Inverted Pendulum 5.1 Material and Methods for 2DOF Experimentations As described in our 2DOF model as in Figure 4.2, two degree of freedom inverted pendulum experimental results also showed two modes of vibration. Figure 5.1 shows the snapshotsofourexperimentalsetupwhilethe2DOFpendulumwasstabilizedinaswinging motion around its upright vertical position under periodic forcing. On the left graph of Figure 5.1 the 2DOF inverted pendulum was under mode 1 vibration; on the right graph of Figure 5.1 the 2DOF inverted pendulum was under mode 2 vibration. Similar to oursDOF experimental setup, our 2DOF apparatus has its base mounted to a stand which is rmly mounted on the center cone of the speaker. In Figure 5.1, the base of the pendulum is labelled as Joint 1, in which the attachment is a ball-bearing made of stainless steel. The 116 rst pendulum rod labelled as L 1 , which is the same sDOF pendulum rod as in Figure 3.1. The rst pendulum rod L 1 measured as L 1 = 80:00:1mm. On top of the rst pendulum rodL 1 weaddedanotherpendulumrodL 2 usinganotherball-bearingjointlabelledasJoint 2tomake 2DOFsetup. ThesecondpendulumrodL 2 measuredasL 2 = 40:00:1mm. All ball-bearing used in our experimental setup are free to rotate on the pendulum XY plane. Each stainless steel made ball-bearing has heavier weight compare with pendulum rods in ordertoful llourtheoreticalmodelasdescribedinFigure4.2. Properlubricationwasused on all ball-bearings before we ran each experiment to ensure that the damping e¤ect due to pendulum joint friction was minimized. We used the same speaker, function generator and ampli ersetupasinsDOFexperimentationsasshowninFigure3.3toour2DOFapparatus. During experiment the speaker cone was in a sinusoidal vertical motion u(t) = "cos!t to move the pendulum base where its forcing amplitude " and forcing frequency ! can be controlled by the function generator and the ampli er. Similar to sDOF experiments, we applied a systematic way to drive the 2DOF pendulum base to appropriate forcing amplitude " and forcing frequency ! while released the pendulum from some initial conditions. All pendulums swinging motion on the XY plane was captured using high speed camera with frame rate 500 frames per second. A ruler with 10mm space tick marks on it was located at a xed location on the right side of the background which we used to calibrate the pendulum physical size from the images capturedfromhighspeedcamera. Wemountedthisrulercarefullytoensurethattherulers straight line was aligned to the vertical direction as close as possible. For convenience, we 117 Figure 5.1: Two degree of freedom inverted pendulum experimental setup. (a)Mode1 vi- bration. (b)Mode2 vibration. 118 rotated the high speed camera with 90 angle which the imageswide aspect ratio could t inmorependulumsswingingmotionincameraviewevenwhenthependulumwasswinging in large angle motion. With new image aspect ratio on 2DOF setup, single image takes 860604 pixels of image resolution which is enough for a 1G byte of computer memory to capture up to 1:7 seconds of data. After we captured a sequence of images taking from high speed camera we used the same image processing tools, NI Vision Builder, as insDOF case to measure 2DOF inverted pendulums motion in a precise way. Figure 5.2 shows a screen shot of the NI Vision Builder used on a PC. Instead of vertical motion as in real experiments,the2DOFinvertedpendulumgraphinFigure5.2showshorizontalmotiondue to 90 rotation of the camera where the gravitational acceleration g is now on a horizontal position pointing to the right. There were many regions of interest shown in Figure 5.2 in a red-color or green-color boxes which indicate the pendulum motion was measured from following the motion of the pendulum shape in the sequence of images captured from high speed camera. Individual pendulum angles can be measured from comparing the the di¤erence of the pendulums straight edge to the straight edge of the ruler mounted on the background. The pendulums response can then be obtained from the evolution of pendulum angles change in di¤erent frame sequences of images. This pendulum response measurements were in high precision since the high speed camera was capturing in a high frame rate as 500 frames per second which resulted in high resolution time traces. 119 Figure 5.2: NI Vision Builder on two degree of freedom experimental setup. 120 5.2 2DOF System Experimental Results Our 2DOF experimentations also showed two modes of vibration. Figure 5.3, 5.4 and 5.5 present our experimental results of 2DOF inverted pendulum under mode 1 vibra- tion; Figure 5.6, 5.7 and 5.8 present our experimental results of 2DOF inverted pendulum under mode 2 vibration. Notice that there always exist the upright vertical mode when the pendulum doesnt not swing as in Figure 5.9, 5.10 and 5.11. We didnt discuss the upright vertical mode much since our research interest was in searching for multiple frequency pe- riodic responses of an inverted pendulum under periodic forcing. However, from stability point of view, the upright vertical mode is the asymptotic solution for the two modes of vibration in a stable condition when the system encounter damping e¤ect. Experimental Results on 2DOF Mode 1 Vibration In Figure 5.3, part (a) shows the input periodic forcing which we present the raw data on top graph and the ltered data at the bottom. We applied the same Butterworth low pass lter as described in sDOF case. The periodic forcing amplitude was " = 0:0085 0:0005m while the periodic forcing frequency was ! = 30:00:1Hz. We started from some initialconditionsawayfromuprightverticalpositionandthis2DOFpendulumwasswinging under mode 1 vibration around its inverted position. In part (b) graph, we also measured the horizontal motion of the pendulum base in order to make sure that there is little e¤ect from the horizontal direction to ful ll our assumption. The horizontal motion measured as y = 0:00:2mm which is negligible compare the amplitude of vertical forcing. In part (c) graph, two inverted pendulum rodsresponse present in time traces with comparison of the 121 input periodic forcing. This shows that two pendulum rods are in mode 1 vibration where their motions are in the same directions. Notice that the damping e¤ect can also be seen from graph in part (c) where two pendulums maximum swinging amplitude decreased in later timing. Graph in part (d) shows the angular velocity of two pendulum rods, their motions also show in a synchronize directions and a detail zoom-in view can be seen from part (f) graph. Part (e) graph shows a zoom-in view of part (c) which indicates that fast time response of two pendulum rods synchronize with the input periodic forcing frequency at ! = 30:00:1Hz. Figure 5.4 presents the same 2DOF experimental results as in Figure 5.3 in a di¤erent manner. We separated Figure 5.4 in left and right portions where left graphs present the response of the rst pendulum rod while right graphs present the response of the second pendulum rod. All graphs were shown in comparison of its raw data points(top) to ltered data(bottom). Notice that in part (e) and (f), the phase portrait generated from raw data points were not recognizable without data ltering. We measured this 2DOF inverted pendulum response as 36 : 1 frequency ratio. Mode 1 pendulum response is very much alike the sDOF case. This can be seen from adding sti¤ness to joint 2 of our 2DOF model to become a sDOF model. The power spectral density graphs of 2DOF inverted pendulum in Figure 5.5 also show similar results as in sDOF case. Top portion of graphs indicate raw experimental data while the bottom graphs are ltered data. There are two clear peaks in our PSD graph in Figure 5.5, the largest peak happened at low frequency around 1Hz which corresponds to the system slow time, the second peak located at around 30Hz which corresponds to periodic forcing frequency. 122 Figure 5.3: Experimental result on 2DOF inverted pendulum under periodic forcing ampli- tude " = 0:0085 0:0005m and periodic forcing frequency ! = 30:0 0:1Hz. System is undermode1vibration. (a)periodicforcing(b)pendulumbasehorizontalmovement(c)time traces of system response (d)time traces of system angular velocity (e)time traces zoom-in view of system response (f)time traces zoom-in view of system angular velocity. 123 Figure 5.4: Experimental result on 2DOF inverted pendulum under periodic forcing ampli- tude " = 0:0085 0:0005m and periodic forcing frequency ! = 30:0 0:1Hz. System is under mode 1 vibration. Left graphs indicate the system response of the rst pendulum; right graphs indicate the system response of the second pendulum. 124 Figure 5.5: Power spectral density graphs of 2DOF inverted pendulum under periodic forcing amplitude " = 0:00850:0005m and periodic forcing frequency ! = 30:00:1Hz. System is under mode 1 vibration. (a)PSD for the rst pendulum (b)PSD for the second pendulum. Experimental Results on 2DOF Mode 2 Vibration In Figure 5.6, part (a) shows the input periodic forcing which we present the raw data on top graph and the ltered data at the bottom. We applied the same Butterworth low pass lter as described in sDOF case. Notice that our input forcing shows periodic pattern in part (a) graph, however, it is not a harmonic function as u(t) = "cos!t. There were two forcing amplitudes appeared in the time traces. This phenomenon happened in all our 2DOF mode 2 experimentations. One possible explanation on this phenomenon is that it may exceed our speakers driving capability while the 2DOF inverted pendulum is under mode 2 vibration. Due to this e¤ect we enlarge the measuring uncertainty on forcing amplitude as " = 0:0080:003m with the periodic forcing frequency as ! = 30:00:1Hz. In part (b) graph, we also measured the horizontal motion of the pendulum base in order to make sure that there is little e¤ect from the horizontal direction to ful ll our assumption. 125 The horizontal motion measured as y = 0:0 0:2mm which is negligible compare the amplitude of vertical forcing. In part (c) graph, two inverted pendulum rodsresponse present in time traces with comparison of the input periodic forcing. This shows that two pendulum rods are in mode 2 vibration where their motions are in the opposite directions. Unlike mode 1 case which has obvious damping e¤ect, mode 2 vibration only show little damping e¤ect as can be seen from part (c) graph. Graph in part (d) shows the angular velocity of two pendulum rods, their motions also show in opposite directions and a detail zoom-in view can be seen from part (f) graph. Part (e) graph shows a zoom-in view of part (c). Figure 5.7 presents the same 2DOF experimental results as in Figure 5.6 in a di¤erent manner. We separated Figure 5.7 in left and right portions where left graphs present the response of the rst pendulum rod while right graphs present the response of the second pendulum rod. All graphs were shown in comparison of its raw data points(top) to ltered data(bottom). Notice that in part (e), the phase portrait generated from raw datapointswerenotrecognizablewithoutdata ltering. Ourexperimentalresultsin 2DOF mode 2 vibration case doesnt match with our numerical simulation results due to di¤erent forcing inputs. Under this usual forcing input, our experimental result showed frequency ratio near 2 : 1 which we couldnt nd from any of our simulation result provided only harmonic forcing input. The power spectral density graphs of 2DOF inverted pendulum in Figure 5.8 show multiple peaks in various frequency locations. Top portion of graphs indicate raw experi- mental data while the bottom graphs are ltered data. 126 Figure5.6: Experimentalresulton2DOFinvertedpendulumunderperiodicverticalforcing. The forcing frequency was ! = 28:00:1Hz. System is under mode 2 vibration. Notice that the input forcing was not harmonic which had two forcing amplitudes involved. 127 Figure5.7: Experimentalresulton2DOFinvertedpendulumunderperiodicverticalforcing. Theforcingfrequencywas! = 28:00:1Hz. Systemisundermode2vibration. Leftgraphs indicatethesystemresponseofthe rstpendulum;rightgraphsindicatethesystemresponse of the second pendulum. 128 Figure 5.8: Power spectral density graphs of 2DOF inverted pendulum under periodic vertical forcing. The forcing frequency was ! = 28:0 0:1Hz. System is under mode 2 vibration. (a)PSD for the rst pendulum (b)PSD for the second pendulum. Experimental Results on 2DOF Upright Vertical Mode Figure 5.9 presents the 2DOF pendulum stabilized in its upright vertical position. This can be obtained from assigning the initial conditions as close to zero as possible at start up. In part (a) shows the input periodic forcing which we present the raw data on top graph and the ltered data at the bottom. We applied the same Butterworth low pass lter as described in sDOF case. In part (b) graph, we measured the horizontal motion of the pendulum base as y = 0:00:1mm. In part (c) graph, two inverted pendulum rods response present in time traces with comparison of the input periodic forcing. This graph indicates that two pendulum rods were in near vertical position without swinging. Same e¤ect can be seen from part (d) where two pendulum rodsangular velocity were small. Same conclusions can be seen in zoom-in view in part (e) and (f). Figure 5.10 presents the same 2DOF experimental results as in Figure 5.9 in a 129 Figure5.9: Experimentalresulton2DOFinvertedpendulumunderperiodicverticalforcing with 0 initial condition. The forcing frequency was! = 30:00:1Hz. Notice that the input forcing was not harmonic which had two forcing amplitudes involved. 130 di¤erent manner. We separated Figure 5.10 in left and right portions where left graphs present the response of the rst pendulum rod while right graphs present the response of the second pendulum rod. All graphs were shown in comparison of its raw data points(top) to ltereddata(bottom). Noticethatinpart(e)and(f),theirphaseportraitgraphsindicate that two pendulum rods stayed in their upright vertical position with very little movement. Underuprightverticalmodeof2DOFinvertedpendulum,itspowerspectraldensitygraphs as in Figure 5.11 cannot conclude any harmonic peak. 131 Figure 5.10: Experimental result on 2DOF inverted pendulum under periodic vertical forc- ing with 0 initial condition. The forcing frequency was ! = 30:00:1Hz. Notice that the input forcing was not harmonic which had two forcing amplitudes involved. 132 Figure 5.11: Power spectral density graphs of 2DOF inverted pendulum under periodic vertical forcing with 0 initial condition. The forcing frequency was ! = 30:00:1Hz. 133 Chapter 6 Summary and Future Work Our research shows that under stable conditions near a pendulums inverted state, as the number of degree of freedom increase, the maximum angle of the pendulums swing- ing motion decreases. In this paper, we presented four ways to characterize the pendulums periodic responses, allowing us to conclude that when the periodic forcing amplitude in- creases,thefrequencyratiodecreasesandthesystemslowtimeresponsewillmovetohigher frequency resulting in faster system response. The power spectral density graph has fewer harmonic peaks when the periodic forcing amplitude is low. As we increase the periodic forcingamplitude, thosefewharmonicpeakswillsplitintomultiplefrequencypeaks. These statements are true for both sDOF case and 2DOF case. In searching for the periodic solutions of an inverted pendulum under periodic vertical forcing to its base, our research shows that our numerical simulation results match reasonably well with our experimental results, quantitatively for the single degree of free- dom case, and qualitatively for the multi-degree of freedom case. Although we carefully 134 constructed our 2DOF experimental setup in order to ful ll our theoretical model as de- scribed in Figure 4.2, each pendulums center of mass still has some o¤set location away from the pivot point which contradicts our initial assumption of point-mass pendulum model. This can be seen from Figure 5.1, in which we measured the combined center of mass for the rst pendulum located at R 1 = 53:7 0:1mm with combined weight mea- sured as m 1 = 3:80:1g; the combined center of mass for the second pendulum located at R 2 = 25:10:1mm with combined weight measured as m 2 = 1:90:1g. We can construct a modi ed theoretical model as described in Figure 6.1 to further match our 2DOF experi- mental setups. We leave this to our future work to tighten up the agreement between our theoretical predictions and experiments. Figure 6.1: Modi ed 2DOF model has each pendulums center of mass o¤set to a new location di¤erent than the pivot point. 135 Bibliography [1] Acheson, D.J., Multiple-nodding oscillations of a driven inverted pendulum. Proc. R. Soc Lond. A448, 89-95. 1995. [2] Acheson, D.J., A pendulum theorem. Proc. R. Soc Lond. A443, 239-245. 1993. [3] Acheson, D.J.&Mullin, T., Upside-down pendulum.Nature, Lond.366, 215-216.1993. [4] Acheson, D.J. & Mullin, T., Ropy magic. New Scientist. February, 32-33. 1998. [5] Arnold,V. I.,Koz lov, V. V., and Neishtadt, A. I., Mathematical Aspects of Classical and Celestial Mechanics, volume 3 of Encyclopaedia of Mathematical Sciences(Berlin: Springer-Verlag). 1985. [6] Arnscott,F.M., Latent roots of tri-diagonal matrices. Proc. Edinburgh Math. Soc., 12. Edinburgh Math. Notes No. 44. 1961. [7] Baillieul,J.,Stableaveragemotionsofmechanicalsystemssubjecttoperiodicforcing.In DynamicsandControlofMechanicalSystems: TheFallingCatandRelatedProblems: Fields Institute Communications (Providence, RI: American Mathematical Society), pp. 1~23. 1993. [8] Baillieul,J., Energy methods for the stability of bilinear systems with oscillatory inputs. Int. J. Robust and Nonlin. Cont., July, Special Issue on the Control of Nonlinear Mechanical Systems, pp. 285~301. 1995. [9] Baillieul, J., The geometry of controlled mechanical systems. In J. Baillieul and J. C. 1998. [10] Baillieul, J., and Weibel, S. P., Scale dependence in the oscillatory control of micro- mechanisms.Submittedtothe1998IEEEConferenceonDecisionandControl(Tampa, FL). 1998. [11] BerliozA,DufourR.&SinhaS.C,ExperimentalandNumericalInvestigationofaNon- linear Autoparametric System. ASME Design Technical Conferences, DETC97/VIB- 4027. 1997. [12] Blackburn, J. A., Smith, H. J. T. & Gronbech-Jensen, N., Stability and Hopf bifurca- tions in an inverted pendulum. Am. J. Phys. 60, 903-908. 1992. 136 [13] Bogoliubov, N. N., and Mitropolsky, Y. A., Asymptotic Methods in the Theory of Non- linear Oscillators. International Monographs on Advanced Mathematics and Physics (New York: Gordon and Breach). 1961. [14] ByrdP.&FriedmanM., Handbook of Elliptic Integrals for Engineerings and Physicists. Springer Verlag. 1954. [15] Champneys, A. & fraser, W. The Indian rope trickfor a continuously exible rod; linearized analysis. Proc. Roy. soc. Lond. A 456, 553-570. 2000. [16] Dombovari Z., Wilson R.E., Stepan G., Estimates of the Bistable Region in Metal Cut- ting. (1)Department of Applied Mechanics, Budapest university of Technology and Economics, Budapest, Hungary. (2)Department of Engineering Mathematics, Univer- sity of Bristol, Bristol, UK. [17] Feng,Z.C.,andLeal,L.G.,Nonlinearbubbledynamics.Ann.Rev.FluidMechanics,29, 201~243. 1997. [18] Flashner, H. & Hsu, C.S., A study of nonlinear periodic systems via the point mapping method. International Journal for Numerical Methods in Engineering. Vol. 19. 185-215. 1983. [19] Gani, Ferawati., Parametric excitation of stay-cables. University of Ottawa (Canada) ; AAT MR01476. 2004. [20] Henry Dan, Perturbation of the Boundary in Boundary-Value Problems of Partial Dif- ferential Equations. Cambridge University Press. ISBN-10 0-521-57491-9. 2005. [21] HowellG.W,OnoK., The connection between n-link chains and the continuum models. Boston University Center for Biodynamics. 1998. [22] Howell G.W, Analysis and control of Super-Articulated Biped Robots. PhD thesis, Boston University. 2000. [23] Je¤rey, A., Handbook of Mathematical Formulas and Integrals. Academic Press Inc. ISBN 0-12-382580-6, 231-234. 1995. [24] Johnson, R.S., Singular Perturbation Theory. Springer Science + Business Media, Inc. ISBN 0-387-23200-1. 2004. [25] Kalmus, H. P., The inverted pendulum. Am. J. Phys. 38, 874-878. 1970. [26] Kaajakari,VillePekka, Ultrasonic surface micromachine actuation: Applications to re- lease, microstructure assembly, and micromotors.The University of Wisconsin - Madi- son ; AAT 3060427. 2002. [27] KhalilHassanK.NonlinearSystems.ThirdEdition.PrenticeHall.ISBN0-13-067389-7, 2002. 137 [28] Kim, S. and Hu, B. Bifurcations and transitions to chaos in an inverted pendulum. Physical Review E, Vol. 58, No. 3, pp. 3028-3035, 1998. [29] Kuo, A.D., Donelan, J.M., and Ruina, A. Energetic consequences of walking like an inverted pendulum: step-to-step transitions. Exerc. Sport Sci. Rev., Vol. 33, No. 2, pp. 88-97, 2005. [30] Kuznetsov Yuri A. Elements of Applied Bifurcation Theory. Third Edition. Springer- Verlag New York, LLC. ISBN 0-387-21906-4, 2004. [31] Levi,M., andWeckesser,W., Stabilization of the inverted pendulumby high frequency vi- brations. SIAM Review, 27, 219~223. 1995. [32] Marcus,M., and Minc,H., A Survey of Matrix Theory and Matrix Inequalities (New York: Dover). 1994. [33] McLaughlin, J. Period-doubling bifurcations in the parametrically forced pendulum. Journal of Statistical Physics, Vol. 24, pp. 375-338, 1981. [34] Meirovitch L, Fundamentals of Vibrations. McGraw-Hill Companies. ISBN 0-07- 288180-1. 2001. [35] MennemR.C., Parametrically excited vibrations in spiral bevel geared systems.Purdue University ; AAT 3150802. 2004. [36] Michaelis M., Stroboscopic Study of the Inverted Pendulum. Am. J. Phys. 53(11) No- vember 1985. [37] Mullin T., Champneys A., Barrie Fraser W., Galan J., Acheson D.. The Indian rod trickvia parametric excitation. Proc. Roy. soc. Lond. [38] Morrison T., Three Problems in Nonlinear Dynamics with 2:1 Parametric Excitation., Cornell University, 146 pages; AAT 3227288. 2006. [39] Nayfeh, A. H., Introduction to Perturbation Techniques. John Wiley & Son, Inc. ISBN 0-471-31013-1. 1981. [40] Nayfeh, A. H., Perturbation Mathods. John Wiley & Son, Inc. ISBN 0-471-63059-4. 1973. [41] Nayfeh, A. H., Problems in Perturbation. John Wiley & Son, Inc. ISBN 0-471-82292-2. 1985. [42] Newton Paul K., The N-Vortex problem: analytical techniques. Applied mathematical sciences. Springer-Verlag New York, Inc. ISBN 0-387-95226-8. 2001. [43] Ng, Leslie Ho-Ming, Three problems in parametric excitation. ,Cornell University, 159 pages; AAT 3050435. 2002. [44] Pai Frank. Highly Flexible Structures: Modeling, Computation, and Experimentation. AmericanInstituteofAeronauticsandAstronautics,Inc.ISBN-101-56347-917-6.2007. 138 [45] Rand R.H, Armbruster D., Perturbation methods, bifurcation theory and computer algebra. Applied Mathematical Sciences, 65. Springer-Verlag, NY, USA. 1987. [46] Rand R. H., Morrison T.M., 2:1:1 Resonance in the Quasiperiodic Mathieu Equation. Nonlinear Dynamics 40:195-203. 2005. [47] Rand R. H., Barcilon A., Morrison T. M. Parametric Resonance of Hopf Bifurcation, Nonlinear Dynamics 39:411-421. 2005. [48] Seredynska M., Hanyga A., Nonlinear di¤erential equations with fractional damping with applications to the 1dof and 2dof pendulum. Acta Mechanica. Vol.176, Iss. 3-4; p. 169. 2005. [49] Smith H. J. T. & Blackburn J. A, Experimental study of an inverted pendulum. Am. J. Phys. 60, 909-911. 1992. [50] Spears, Brian Keith, Nonlinear dynamics of quasiperiodically driven oscillators: Knot- ted torus attractors and their bifurcations.UniversityofCalifornia,Berkeley,115pages; AAT 3147013. 2004. [51] Stephan G., Insperger T. & Szalai R., Nonlinear Dynamics of High-speed Milling, Budapest University of Technology and Economics, Hungary. 2005. [52] StephanG.,InspergerT.&SzalaiR., Delay, Parametric Excitation, and the Nonlinear DynamicsofCuttingProcess.InternationalJournalofBifurcationandChaos159:2783. 2005. [53] StephensonA.,MemoirsandProceedingsoftheManchesterLiteraryandPhilosophical Society 52(8), 1-10. On a new type of dynamic stability. 1908. [54] Stephenson A., Philosophical Magazine 15, 233-236. On induced stability. 1908. [55] Strogatz, S.H., Nonlinear Dynamics and Chaos. Addison-Wesley Publishing Company. ISBN 0-201-54344-3. 1994. [56] Weibel S., Baillieul J., Lehman B., Equilibria and stability of an n-pendulum forced by rapid oscillations, Decision and Control, Proceedings of the 36th IEEE Conference on Volume2, Issue , 10-12 Dec 1997 Page(s):1147-1152 vol.2, Digital Object Identi er 10.1109/CDC.1997.657602. 1997. [57] Weibel S., Applications of Qualitative methods in the nonlinear control of superarticu- lated mechanical systems. PhD thesis, Boston University. 1997. [58] Weibel S., Baillieul J., Oscillatory control of bifurcations in rotating chains. In Pro- ceedings of 1997 ACC, Albuquerque, NM, pp. 2713~2717. 1997. [59] Weibel S., Baillieul J., Averaging and energy methods for robust open-loop control of mechanical systems. In J. Baillieul, S. Sastry, and H. Sussmann (Eds), Essays in Math- ematical Robotics, (Berlin: Springer-Verlag), pp. 203~270. 1998. 139 [60] Weibel S., Baillieul J., Kaper T. J., Small-amplitude periodic motions of rapidly forced mechanical systems.InProceedingsof34thIEEECDC,NewOrleans,LA,pp.533~539. 1995. [61] Weibel,S.,Baillieul,J.,Lehman,B., Equilibria and stability of an n-pendulum forced by rapid oscillations. In Proceedings of 36th IEEE CDC, San Diego, CA, pp.1147~1152. 1997. [62] Weibel S., Baillieul J., Open-loop oscillatory stabilization of an n -pendulum. Interna- tional Journal of Control, Volume 71, Number 5, 16 November 1998 , pp.931-957(27). [63] Weibel S.,Kaper T. J., Baillieul J., Global dynamics of a rapidly forced cart and pen- dulum. Nonlin. Dyn., 13, 131~170. 1997. [64] P. A. Vela, Averaging and Control of Nonlinear Systems. Ph.D. thesis, California In- stitute of Technology, 2003. [65] P. A. Vela, K. A. Morgansen, and J. W. Burdick. Second-order averaging methods for oscillatory control of underactuated mechanical systems. In Proceedings of the Ameri- can Control Conference, 2002. [66] Yagasaki K., Ichikawa T. Higher-order averaging for periodically forced weakly nonlin- ear systems. International Journal of Bifurcation and Chaos, 9(3):51931, 1999. [67] Zounes R. S, Rand R. H., Subharmonic resonance in the non-linear Mathieu equation, Int. J. Nonlin. Mech. 37, 43. 2002. 140 Appendix Below are power series expansions for transition curves in Mathieus equation cal- culated from solving Hills equation using harmonic balance in Rand [47]. Those transition curves data shown in Figure 2.2, 2.3 and 2.4 were calculated using these power series ex- pansions. (n = 0) = " 2 2 + 7" 4 32 29" 6 144 + 68687" 8 294912 123707" 10 409600 + 8022167579" 12 19110297600 +::: (6.1) (n = 1) 1 l = 1 4 " 2 8 + " 3 32 " 4 384 11" 5 4608 + 49" 6 36864 55" 7 294912 83" 8 552960 + 12121" 9 117964800 114299" 10 6370099200 192151" 11 15288238080 + 83513957" 12 8561413324800 ::: 1r = 1 4 + " 2 8 " 3 32 " 4 384 + 11" 5 4608 + 49" 6 36864 + 55" 7 294912 83" 8 552960 12121" 9 117964800 114299" 10 6370099200 + 192151" 11 15288238080 + 83513957" 12 8561413324800 +::: (6.2) (n = 2) 2 l = 1 " 2 12 + 5" 4 3456 289" 6 4976640 + 21391" 8 7166361600 2499767" 10 14447384985600 + 1046070973" 12 97086427103232000 ::: 2r = 1+ 5" 2 12 763" 4 3456 + 1002401" 6 4976640 1669068401" 8 8=7166361600 + 4363384401463" 10 14447384985600 40755179450909507" 12 97086427103232000 +::: (6.3) 141 (n = 3) 3 l = 9 4 + " 2 16 " 3 32 + 13" 4 5120 + 5" 5 2048 1961" 6 1474560 + 609" 7 3276800 + 4957199" 8 33030144000 872713" 9 8493465600 + 421511" 10 23488102400 + 16738435813" 11 1331775406080000 572669780189" 12 58706834227200000 +::: 3r = 9 4 + " 2 16 + " 3 32 + 13" 4 5120 5" 5 2048 1961" 6 1474560 609" 7 3276800 + 4957199" 8 33030144000 + 872713" 9 8493465600 + 421511" 10 23488102400 16738435813" 11 1331775406080000 572669780189" 12 58706834227200000 +::: (6.4) (n = 4) 4 l = 4+ " 2 30 + 433" 4 216000 5701" 6 170100000 112236997" 8 31352832000000 + 8417126443" 10 123451776000000000 + 2887659548698709" 12 265470699110400000000000 +::: 4r = 4+ " 2 30 317" 4 216000 + 10049" 6 170100000 93824197" 8 31352832000000 + 21359366443" 10 123451776000000000 2860119307587541" 12 265470699110400000000000 +::: (6.5)
Abstract (if available)
Abstract
A pendulum is statically unstable in its upright inverted state due to the Earth's gravitional attraction which points downward. However, with proper forcing, the pendulum can be stabilized in its upright inverted state. Special interest is on periodic vertical forcing applied to the pendulum's base to stabilize it around the upright inverted equilibrium. Many researchers have studied how to stabilize the system by varying various parameters, in particular its amplitude and frequency. Most have focused on the single degree of freedom inverted pendulum case, which with linear assumption can be described via Mathieu's equation. The system stability can then be characterized by Floquet theory. Our focus is on searching for the periodic solutions inside the linearly stable region of the pendulum's inverted state when the pendulum is under proper periodic forcing. Our research shows that under appropriate excitation by controlling the forcing amplitude and frequency, the pendulum can maintain certain periodic orbits around its inverted state which we characterize in a systematic way.
Linked assets
University of Southern California Dissertations and Theses
Conceptually similar
PDF
New approaches in modeling and control of dynamical systems
PDF
Modeling and analysis of parallel and spatially-evolving wall-bounded shear flows
PDF
Computationally efficient design of optimal strategies for passive and semiactive damping devices in smart structures
PDF
Modeling and dynamic analysis of coupled structure-moving subsystem problem
Asset Metadata
Creator
Chen, Cheng-Yuan Jerry
(author)
Core Title
Multiple degree of freedom inverted pendulum dynamics: modeling, computation, and experimentation
School
Viterbi School of Engineering
Degree
Doctor of Philosophy
Degree Program
Aerospace
Publication Date
05/18/2009
Defense Date
12/10/2008
Publisher
University of Southern California
(original),
University of Southern California. Libraries
(digital)
Tag
floquet theory,inverted pendulum,Mathieu's equation,nonlinear dynamics,OAI-PMH Harvest,parametric excitation
Language
English
Contributor
Electronically uploaded by the author
(provenance)
Advisor
Newton, Paul K. (
committee chair
), Flashner, Henryk (
committee member
), Jonckheere, Edmond A. (
committee member
), Redekopp, Larry G. (
committee member
)
Creator Email
chenchen@usc.edu,chenchen1112@sbcglobal.net
Permanent Link (DOI)
https://doi.org/10.25549/usctheses-m2262
Unique identifier
UC1284983
Identifier
etd-Chen-2707 (filename),usctheses-m40 (legacy collection record id),usctheses-c127-238695 (legacy record id),usctheses-m2262 (legacy record id)
Legacy Identifier
etd-Chen-2707.pdf
Dmrecord
238695
Document Type
Dissertation
Rights
Chen, Cheng-Yuan Jerry
Type
texts
Source
University of Southern California
(contributing entity),
University of Southern California Dissertations and Theses
(collection)
Repository Name
Libraries, University of Southern California
Repository Location
Los Angeles, California
Repository Email
cisadmin@lib.usc.edu
Tags
floquet theory
inverted pendulum
Mathieu's equation
nonlinear dynamics
parametric excitation