DESIGN AND TUNING OF CENTRIFUGAL PENDULUM VIBRATION ABSORBERS FOR NONLINEAR RESPONSE By Mustafa Ali Acar A DISSERTATION Submitted to Michigan State University in partial fulfillment of the requirements for the degree of Mechanical Engineering – Doctor of Philosophy 2017 ABSTRACT DESIGN AND TUNING OF CENTRIFUGAL PENDULUM VIBRATION ABSORBERS FOR NONLINEAR RESPONSE By Mustafa Ali Acar Centrifugal pendulum vibration absorbers (CPVAs) are used to reduce engine-order torsional vibrations in rotating systems. They are widely used in light aircraft piston engines and helicopter rotors and have recently been introduced for smoothing torsional vibrations in automotive powertrain applications. These absorbers make use of the centrifugal field due to rotation, in place of elastic elements, for their restoring force, so that they are tuned to a particular engine order, rather than a particular frequency, making them ideally suited for these applications. Increasing demands on fuel economy have led to engine downsizing and downspeeding, resulting in harsh torsional excitations being exerted on powertrain components, which are ultimately felt by passengers. In order to maintain durability and NVH performance specifications for these engines, sophisticated vibration control solutions are required, many of which currently involve CPVAs. The contributions made by this study are in three major topics, all related to the design and performance of CPVAs. The first is the development of a new numerical tool for the rapid analysis of the response of complex CPVA systems. This algorithm recasts the equations of motion in a way that significantly extends the applicability and accuracy of the widely known harmonic balance method. This method essentially automates the analysis of steady state responses of systems with multiple CPVAs, including their stability characteristics. This method is much faster than brute-force numerical simulations, and it allows designers to explore the response of systems that are not amenable to more traditional analysis tools, such as perturbation methods. The capabilities of this new tool are demonstrated through two example systems that are known to exhibit rich dynamical characteristics. The second topic is the investigation of a new configuration that uses CPVAs to eliminate crankshaft torsional resonances under order excitation. This dynamical system is a combination of a frequency based element, a shaft torsional vibration mode, and an order based element oscillator, the CPVA. We show that the linear natural frequencies of this system undergo eigenvalue veering as the engine speed is varied near the resonance point. The structure of this veering suggests that with proper tuning of the absorber, one can eliminate the shaft torsional resonance. We use perturbation methods to show that these results extend to operating conditions where the CPVA amplitudes are large and its response becomes nonlinear. The third topic deals with the design and analysis of a new type of kinematic suspension that increases the effective inertia of the absorber, thereby reducing the packaging space for these absorbers. We show that by employing non-symmetric cutouts for the rollers used in standard bifilar (two point) suspensions of CPVAs, we can specify both the rotational and translation motions of the absorber relative to the rotor. This allows designers to increase the effective inertia of the CPVA, thereby providing better vibration correction for a given amount of absorber mass. The dynamical response characteristics for this system are studied using both perturbation methods and the newly developed harmonic balance method. It is shown that this so-called “rocking” absorber configuration provides an improvement of about 15% when compared to its traditional counterpart. ACKNOWLEDGMENTS I would like to express my gratitude to both of my advisors Professor Steven Shaw and Professor Brian Feeny. It has been an extremely rewarding experience to be a part of their research and to be able work with them. Steve has always been an incredibly supportive, considerate and resourceful advisor and a great friend. Apart from his invaluable contribution to my technical development, I had the first hand experience of observing his amazing way of navigating one through the journey of PhD. Acknowledgment sections of theses of all students Brian had advised would have an expression that he is one the nicest persons one can encounter as an advisor and a friend. While there is absolutely no doubt on that, the questions Brian came up with throughout my research unveiled many possible problems that I could have overlooked. I am grateful for the support I received from Prof. Alan Haddow in every step of the upkeep, update and operation of our experimental setup. I would also like to thank my committee members Professor Ranjan Mukherjee and Professor Hassan Khalil for their support and contributions on my dissertation. I am thankful to Bruce Geist for his continuous support of our research and for availing his unique mathematician-from-the-industry perspective to us. I am very glad to be able keep working with him after graduation. I thank my fellow lab mates and colleagues Scott, Ming, Pavel, Ori, Johannes, Venkat, Abhisek, Xing, Smruti, Rickey, Ayse and Fatemeh for the stimulating discussions we had and the pleasant time we shared over the course of my PhD. I also thank Michael and Jessica for their help on my experimental studies. Special thanks go to my wonderful wife Gizem, who is both comforting and exciting. She is someone I learned from a lot and someone I taught a lot. Few people have the chance of iv being married to one of the brightest minds of their field of study like I had. Beyond question, I benefited from her support and presence throughout the course of graduate school. I am already looking forward to the new adventures we will take together. Last but not least, I would like to thank to my amazing family; my mother Aynur, my father Ziya, my sisters Busra and Tugba. Most of the time we could only see each other virtually, but the support and love I received from them never felt anything less than real. This work was supported by the US National Science Foundation under grant CMMI1100260 as well as by Fiat Chrysler Automobiles (FCA) and by Ford Motor Company. v TABLE OF CONTENTS LIST OF TABLES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viii LIST OF FIGURES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ix CHAPTER 1 INTRODUCTION . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1 The Absorber Path - Background . . . . . . . . . . . . . . . . . . . . . . . . 1.1.1 Tuning of point mass CPVAs . . . . . . . . . . . . . . . . . . . . . . 1 6 9 CHAPTER 2 APPLICATION OF THE HARMONIC BALANCE METHOD TO THE EQUATIONS OF MOTION OF SYSTEMS WITH CPVAS . . . 2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.2 Dynamical Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.2.1 Non-dimensionalization . . . . . . . . . . . . . . . . . . . . . . . . . . 2.2.2 Reforming the equations for HBM . . . . . . . . . . . . . . . . . . . . 2.3 Case Studies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.3.1 Response of Circular Path CPVA with High Inertia Ratio . . . . . . . 2.3.2 Bifurcation to Nonunison Response of Two Identical Tautochronic CPVAs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . CHAPTER 3 LINEAR AND NONLINEAR DYNAMICS OF ON FLEXIBLE ROTOR SYSTEMS . . . . . . 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . 3.2 Dynamical Model . . . . . . . . . . . . . . . . . . . . . 3.2.1 Non-dimensionalization . . . . . . . . . . . . . . 3.3 System Response . . . . . . . . . . . . . . . . . . . . . 3.3.1 Scaling and Averaging . . . . . . . . . . . . . . 3.4 Case Studies . . . . . . . . . . . . . . . . . . . . . . . . 3.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . 30 35 SYSTEMS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 37 40 41 43 43 46 52 CHAPTER 4 INCREASING THE EFFECTIVE INERTIA OF BIFILAR SUSPENSION CPVAS BY EMPLOYING THEIR ROTATIONAL INERTIA THROUGH KINEMATIC ALTERATIONS . . . . . . . . . . 4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2 Cutout Kinematics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3 Dynamical Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.1 Non-dimensionalization . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.2 Absorber Tuning . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.2.1 Linear Tuning . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.2.2 Tautochronic Tuning . . . . . . . . . . . . . . . . . . . . . . 54 54 57 60 62 63 67 68 vi CPVA . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 10 14 16 17 21 22 4.4 4.5 4.3.3 Scaling and Averaging . . . . . . . . . . . . . . . . . . . . . . . . . . Case Studies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 71 77 CHAPTER 5 CONCLUSIONS . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79 BIBLIOGRAPHY . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83 vii LIST OF TABLES Table 2.1: Definition of dimensionless variables. . . . . . . . . . . . . . . . . . . . . . 16 Table 3.1: Definition of the non-dimensional variables. . . . . . . . . . . . . . . . . . 42 Table 4.1: Definition of dimensionless variables. . . . . . . . . . . . . . . . . . . . . . 62 Table 4.2: Parameter sets for the rocking rate study. . . . . . . . . . . . . . . . . . . 73 viii LIST OF FIGURES Figure 1.1: Left: CPVAs installed on an helicopter rotor; from Steve Shaw. Right: CPVAs on a light aircraft engine crankshaft; from [1]. . . . . . . . . . . . 2 Figure 1.2: Illustration of four point mass CPVAs, along with their paths on the rotor [2]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 Figure 1.3: The bifilar suspension of a CPVA used as a test specimen at MSU. . . . 5 Figure 1.4: Circular, tautochronic epicycloidal, and cycloidal paths. The radius of curvature of all path in this formulation start from the same value, which implies that they all have the same linear tuning. Except for the circular path, the radius of curvature gradually decreases as per the path definition, and the location of the center of curvature varies as one moves along the paths. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 Figure 2.1: Schematic of a rotor fitted with CPVAs. . . . . . . . . . . . . . . . . . . 15 Figure 2.2: Harmonic amplitudes of the response of a circular path CPVA system, where N = 1, n ˜ 1 = 2.01, α1 = 0.9, µa,1 = 0.005, H = 48, n = 2. Only the n, 2n, and 3n harmonics, out of 48, are shown. The black curve in (a) indicates the reference rotor response at order n when the absorber is locked. Dashed parts of the curves indicate unstable periodic solutions. 24 Figure 2.3: Maximum Floquet multiplier magnitudes of the system in Figure 2.2 at corresponding periodic solutions. Values above one indicate an unstable solution. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 Figure 2.4: Rotor angle (acting similar to time) domain reconstruction of a periodic solution of the system in Figure 2.2. This solution belongs the lower stable branch at Γ = 0.45. The blue curve is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. (a) Rotor Acceleration, (b) Absorber Displacement . . . . . . . . . . . . . . . . . . 26 ix Figure 2.5: Rotor angle domain reconstruction of a periodic solution of the system in Figure 2.2. This solution belongs the middle unstable branch at Γ = 0.37. The blue curve is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. Note that even numerical solution and integration tolerances are sufficient to cause an ultimate deviation from the periodic solution. In this case, the absorber response grew to a point that the numerical integration crashed. (a) Rotor Acceleration, (b) Absorber Displacement . . . . . . . . . . . . . . . . . . . . . . . . . . 27 Figure 2.6: Rotor angle domain reconstruction of a periodic solution of the system in Figure 2.2. This solution belongs the middle unstable branch at Γ = 0.20. The blue curve is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. In this case, the numerically integrated solution eventually jumps down to the lower stable branch. (a) Rotor Acceleration, (b) Absorber Displacement . . . . . . . . . . . . . . . . . . 28 Figure 2.7: Rotor angle domain reconstruction of a periodic solution of the system in Figure 2.2. This solution belongs the upper stable branch at Γ = 0.15. The blue curve is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. (a) Rotor Acceleration, (b) Absorber Displacement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 Figure 2.8: Harmonic amplitudes of the response of the rotor with two identical tautochronic CPVAs, with N = 2, n ˜ 1,2 = 1.50, α1,2 = 0.1, µa,1,2 = 0.005, H = 16, n = 1.5 and the order of the first assumed harmonic is 0.5. (a) Rotor Acceleration, (b) Absorber Displacement . . . . . . . . . . 32 Figure 2.9: Rotor angle domain reconstruction of a periodic solution of the system in Figure 2.8. The excitation torque Γ = 0.12 in this case correspond to a solution where the two CPVAs move in unison. The blue curve in (a) is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. Numerical integration results of the CPVA signals are indistinguishable from the HBM results and are not shown for the sake of clarity in the figure. (a) Rotor Acceleration, (b) Absorber Displacement . . . . . . . . 33 x Figure 2.10: Rotor angle domain reconstruction of a periodic solution of the system in Figure 2.8. This solution is where the two CPVAs are non-synchronous at Γ = 0.14. The blue curve in (a) is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. As above, numerical integration results of the CPVA signals are not shown. (a) Rotor Acceleration, (b) Absorber Displacement . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 Figure 3.1: A crankshaft with torsional failure [3] . . . . . . . . . . . . . . . . . . . . 38 Figure 3.2: Illustration of the model for a single rotor torsional mode and a single CPVA. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 Figure 3.3: Eigenfrequencies of the linearized system versus rotor speed σ for n ˜ = 2. Blue: α → 0, Green: α = 0.005, Red: α = 0.01 . . . . . . . . . . . . . . . 44 Figure 3.4: Tautochronic path steady state CPVA (a1 ) and rotor deflection (a2 ) steady state amplitude values. Blue: stable, Red: unstable, Red-dashed: ˆ = 0.01 to 0.85. . . . . . . CPVA cusp amplitude. n = 2, = 0.1, φ = 0, Γ 49 Figure 3.5: Softening path steady state CPVA (a1 ) and rotor deflection (a2 ) steady state amplitude values. Blue: stable, Red: unstable, Red-dashed: CPVA ˆ = 0.01 to 0.85. . . . . . . . . . cusp amplitude. n = 2, = 0.1, φ = −5, Γ 50 Figure 3.6: Hardening path steady state CPVA (a1 ) and rotor deflection (a2 ) steady state amplitude values. Blue: stable, Red: unstable, Red-dashed: CPVA ˆ = 0.01 to 0.85. . . . . . . . . . cusp amplitude. n = 2, = 0.1, φ = 5 , Γ 51 Figure 4.1: Illustration of different types of CPVA suspension designs. (a) Ring/Roller Type, (b) Compound Pendulum Type, (c) Bifilar Suspension Type . . . 54 Figure 4.2: Conceptual drawing of the system. Here the black outlined boxes represent the motion of the absorber for a standard bifilar mechanism while the red dashed boxes represent the positions of the rocking absorber. S is the arclength coordinate that indicates the position of the COM of the absorber, and ϕ(S) is the rocking angle of the absorber as a function of S. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55 Figure 4.3: Illustration of moving a mass over a flat surface with a roller. . . . . . . 57 Figure 4.4: Vectors indicating that each point on the absorber will have different velocities as a result of the combined translation of the COM and rotation about the COM. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58 xi Figure 4.5: Samples of cutouts for rocking bifilar absorbers for the same COM path with different rocking rates. Here, the red and green dashed lines represent cutouts on the rotor while red and green solid lines represent the absorber cutouts. The solid outlined circles represent the rollers shown at their center position, and the blue dashed line is the COM path. The rocking rates, moving from top left to bottom right, are set such that the rocking angle at the cusp becomes 0◦ , 7◦ , 14◦ , 21◦ . . . . . . . . . . . 59 Figure 4.6: Rocking and tuning values of the CPVA for the cases in Table 4.2. . . . . 75 Figure 4.7: Harmonic amplitudes of order n of the response of the systems in the cases from Table 4.2. The black curve in (a) indicates the reference rotor response at order n, when the absorber is locked at its vertex. . . . 76 xii CHAPTER 1 INTRODUCTION Internal combustion engines are subject to torsional vibrations caused by the periodic changes in cylinder pressures during the stroke cycle. These pressure changes are reflected as oscillations in the generated torque through the kinematics of the slider-crank mechanism. The rate at which these oscillations occur follows the crank angular speed, since the slider-crank kinematics and the combustion cycle do not change with speed. The resulting excitation is referred to as engine order excitation, that is, excitation that has a fixed number of cycles per revolution, rather than a fixed frequency. These fluctuations on the generated torque are unwanted, probably without exception, and systems using combustion engines need to be designed by taking this into account. In a broad sense, the severity of the generated torsional oscillations and the weight restrictions on the system have competing effects on the outcome of the design. This is why generally Diesel engines are heavier, as the components of them undergo harsher combustion cycles when compared to gasoline engines. In the past, it was mostly small aircraft engines that needed to simultaneously survive under severe order excitation and be lightweight [4]. Modern automotive engines are now required to achieve similar capabilities, due to weight and fuel economy considerations. The energy efficiency of an automobile engine can be significantly improved when designed to operate at higher than atmospheric intake pressures. This practice ensures higher combustion efficiency and allows that a given amount of power can be generated at lower speeds, thus reducing pumping losses [5]. However, this also leads an increase in the torsional vibrations generated by the engine. Therefore increasing demands for downsizing and downspeeding engines necessarily require additional solutions for reducing torsional vibrations. Solutions that are generally employed to deal with said torsional vibrations include single 1 or dual mass flywheels, energy inefficient torque converter locking calibrations, as well as centrifugal pendulum vibration absorbers (CPVAs). Figure 1.1: Left: CPVAs installed on an helicopter rotor; from Steve Shaw. Right: CPVAs on a light aircraft engine crankshaft; from [1]. Reduction of torsional vibrations through a CPVA is a dynamic effect and is fundamentally different from the similar gains that can be achieved by adding more rotational inertia to the rotor, which is undesirable for both weight and responsiveness considerations. The fact that CPVAs require small, and in some cases no extra, mass addition to the rotating system, yet reduce torsional vibration levels for any operating speed as a passive device, render them good candidates for systems with demanding performance requirements [6, 7, 8, 9]. Design parameters that control the behavior of an individual CPVA are basically its inertial parameters and kinematics. Kinematics set the absorber dynamic tuning by dictating how the CPVA mass will move relative to its host structure (see Figure 1.2). Depending on the amplitude of this motion, nonlinear effects can emerge. While for small amplitude motion of the mass, the linearized equations of motion provide sufficiently accurate results for design, the validity of these diminish as the amplitude grows and, more importantly, phenomena that linear analysis cannot predict may be encountered. Thus, tuning a CPVA 2 involves creating a entire path of travel, so that it takes response amplitude into account. θ R(S) ρ0 S Figure 1.2: Illustration of four point mass CPVAs, along with their paths on the rotor [2]. Although, the number of parameters that govern a single CPVA is small, the outcome of the system behavior can be quite complicated depending on the configuration and excitation characteristics of a particular system. In the literature, one can find that a wide variety of problems involving CPVAs have been analyzed both theoretically and experimentally. A thorough chapter in [10] was dedicated to CPVA problems and presents a complete review of the state of the art on the subject at the date of its publication (1968), including examples from several aircraft engine applications, theoretical design charts, and detailed design criteria. DenHartog [11] and Newland [12] are among the early studies that emphasize the importance of the nonlinear nature of CPVA systems. Both early practical experiences and subsequent theoretical studies found that circular path CPVAs have a softening characteristic in their tuning order as the amplitude of oscillation grows. This was causing detrimental, and sometimes disastrous, bifurcations in the system response. A patent by Madden [13] proposed that implementing cycloidal paths for the CPVA kinematics resulted in hardening tuning behavior, thus overcoming this problem. The path formulation in [14] generalized 3 circular and cycloidal paths using a control parameter that allowed the path to vary continuously (in design space). These paths included an amplitude independent isolated tuning order, associated with the tautochronic epicyloid. These types of paths are now in wide use [7, 15, 16, 9, 17, 18, 19, 20]. This study elaborates the analytical and numerical investigation of several systems involving CPVAs. Chapter 2 describes the application of a harmonic balance method for the numerical solution of equations governing the steady-state response of systems with multiple CPVAs of general type. Chapter 3 describes the use of CPVAs for suppressing torsional resonances, and chapter 4 provides results from a study in which CPVAs both translate and rotate relative to their host rotor, thereby increasing their effective inertia. These topics are elaborated on in the following paragraphs. The harmonic balance method (HBM) is a powerful analysis tool for nonlinear vibrating systems, provided that the forms of the nonlinearities of the system translate into a manageable set of algebraic equations. The authors of [21] created a framework that modifies the structure of the equations of motion involving a wide variety of nonlinearities into a quadratic polynomial form, which then can be approximated using HBM with as many assumed harmonics as the problem needs for satisfactory accuracy. In this study, we employ this framework for the analysis of CPVAs. The crucial step of this framework is the recasting of the variables into the required form. It has been shown that the dimensionless equations of motion for point mass CPVAs with general paths fitted to a rigid rotor can be put into the quadratic polynomial form. Two benchmark problems with known dynamical characteristics are investigated and the results show that this approach provides a powerful tool for investigating steady-state responses of these absorber systems, including their stability characteristics. This will be very beneficial for parameter studies and design evaluations of CPVA systems that do not allow for the application of perturbation methods, and/or those that make direct numerical simulations very time consuming. Chapter 3 considers the effects of CPVAs on torsionally flexible rotors. More specifically, 4 it considers the possibility of eliminating shaft torsional resonance by the addition of CPVAs. The most widely used purpose of CPVAs, as previously discussed, is to reduce essentially rigid body vibrations of the host structure, the rotor. This means that when analyzing such systems, the rotor is assumed to be a rigid rotating disc. However, in this study, by using a lumped inertia approach, we allow the rotor model to capture one or more torsional resonances of the shaft. As the operating speed of this flexible rotor is varied under order excitation, the frequency of excitation will vary with the speed, thus creating the possibility of matching the excitation frequency with a shaft resonance at certain speeds, inducing a resonant response. Linear modal analysis is used to show that in fact CPVAs can eliminate such resonances altogether, if properly tuned. Then, using perturbation analysis methods, the extent to which this effect is applicable when the absorber amplitudes become finite is investigated. The studies in [17, 22] describe the most similar intention of using CPVAs, although they consider flexural vibrations of a host structure. They included CPVAs to the model of coupled flexible turbine blades with the purpose of eliminating blade resonances. Figure 1.3: The bifilar suspension of a CPVA used as a test specimen at MSU. In Chapter 4 we propose a kinematic variation to the most common way of implementing CPVAs, the bifilar suspension mechanism, with the intention of increasing the effective absorber inertia, and thereby its capacity for reducing vibration for a given amount of mass. 5 The bifilar mechanism uses two sets of cutouts on the absorber and on the rotor, and a pair of rollers, to guide the motion of the absorber along its prescribed path as shown in Figure 1.3. In this configuration, the absorber does not rotate relative to the rotor, and thus its rotational inertia about its center of mass (COM) dynamically acts as an increased inertia of the rotor, while its mass acts like a point mass. The modification included in this study imposes a certain rotational motion to the absorber, relative to the rotor, by judicious adjustments of the cutouts. This allows the rotational inertia of the absorber to join the effort in counteracting the excitation, eventually getting more correction for the same amount of mass, when the rest of the conditions are the same. Since the rate of rotation imposed by this method is controlled by the cutouts, it becomes another absorber design parameter, like the linear and nonlinear tuning path parameters. 1.1 The Absorber Path - Background The theoretical analysis conducted in the subsequent chapters share certain similar techniques in the modeling and solution procedures. The most common aspect of these is the way the absorber motion paths are formulated. This background section serves the purpose of summarizing the path definition used in this and many previous studies. The dynamic behavior of CPVAs depends on the prescribed paths that constrain the motion of their centers of mass (COM) with respect to the rotating frame in which they move. This constraint means that each CPVA in a system has a single degree of freedom. Therefore, a single generalized coordinate is sufficient to describe the motion of each CPVA. The choice of the coordinate does not affect the result of the dynamics, however, choosing a coordinate that simplifies the analysis is of importance. In this study, an arc-length coordinate along the path (denoted by S) is utilized to describe the paths. The farthest point from the center of rotor is defined as the vertex of this arc and the origin of the path coordinate, S = 0, is taken to be at this vertex. 6 -0.11 Tautochronic Epicycloidal Circular Cyloidal -0.115 -0.12 -0.125 -0.13 -0.135 -0.14 -0.03 -0.02 -0.01 0 0.01 0.02 0.03 Figure 1.4: Circular, tautochronic epicycloidal, and cycloidal paths. The radius of curvature of all path in this formulation start from the same value, which implies that they all have the same linear tuning. Except for the circular path, the radius of curvature gradually decreases as per the path definition, and the location of the center of curvature varies as one moves along the paths. In [14], such a path formulation is derived as a family of paths with two parameters which control the linear and nonlinear tuning of the CPVA. Using these parameters, the radius of curvature of the path (ρ) along the arc-length is specified as follows. ρ(S)2 = ρ20 − λS 2 (1.1) where ρ0 is the vertex radius of curvature of the path and λ is a dimensionless real number that determines the gradual reduction in the radius of curvature as one moves along the path. Note that the paths are assumed to be symmetric about S = 0. For λ = 0 this path is a circle, which has constant radius of curvature everywhere, and at λ = 1, it is a cycloid. Values of λ between 0 and 1 yield epicycloidal paths, one of which is quite special, as described below. As described subsequently, in the equations of motion it is convenient to general functions R(S) and G(S) that depend on the path, where R(S) is the distance from the rotor center 7 of rotation to the absorber COM, and G(S) = R2 (S) − 1 4 d R2 (S) dS 2 (1.2) relates to a projection of the tangent of the path onto the moment arm relative to the center of rotor at S. Note that S, ρ(S), ρ0 , R(S) and G(S) are all physical quantities of length, and in the dynamical analysis procedures in following chapters, these are nondimensionalized by their ratio to the vertex distance of the absorber to the rotor center, i.e. R(0). As the analysis is conducted in terms of nondimensional quantities, it is also convenient to prescribe the path kinematics, which determine the absorber tuning, in nondimensional form. One can always trace back the actual physical dimensions using the scaling factor, R(0). When using the path formulation Equation (1.1), there is an exact closed form expression for R(S) in terms of ρ0 and λ, as well as for its nondimensional counterpart r(s) = R(R(0)s)/R(0)). However, it is convenient to use the following Taylor series expansion of r(s) in the analysis, r2 (s) = 1 + φ0 s2 + φ1 s4 , (1.3) where φ0 = 1 − ρˆ−1 0 φ1 = − (1.4) −1 + λ2 + ρˆ0 12ˆ ρ30 with ρˆ0 := ρ0/R(0). (1.5) Here, φ0 controls small amplitude (linear) behavior of the absorber and φ1 controls nonlinear characteristics. This expansion provides a good approximation for many useful paths in the Denman family Equation (1.1), since s = S/R(0) is always less than unity. Note that, as long as the definition of G(S) (and its dimensionless version g(s)) is kept in its exact form Equation (1.2) in the equations of motion, this polynomial expression for r(s) does not compromise the accuracy of the analysis for these paths. 8 Unlike a simple circular path, these paths have a finite range that is set by a cusp point at which g(scusp ) = 0. For general paths simple numerical root finding algorithms can calculate scusp such that g(scusp ) = 0. For φ1 = 0, this cusp point is easily determined to be [7] scusp ≈ 1.1.1 1 −φ0 (1 − φ0 ) . (1.6) Tuning of point mass CPVAs The most common way of implementing a CPVA design is to employ rollers inside two symmetric, inverted, cutouts on the absorber and the rotor, resulting in the common bifilar suspension mechanism; see Figure 1.3. This configuration provides excellent packaging features and is the most widely used suspension in practice. The standard version of this suspension ensures that the absorber only translates with respect to the rotor and essentially renders it a point mass, in terms of its absorbing inertia. In this case, constructing the absorber path for a desired value of linear tuning n ˜ becomes easy. Specifically, the quadratic coefficient of the absorber path, φ0 , dictates n ˜ , as φ0 = −˜ n2 . (1.7) The nonlinear part of the path is governed by φ1 (or, equivalently, λ), as described in [7, 15, 16, 9, 17, 18, 19, 20], and allows for a variety of nonlinearly softening and hardening paths. 9 CHAPTER 2 APPLICATION OF THE HARMONIC BALANCE METHOD TO THE EQUATIONS OF MOTION OF SYSTEMS WITH CPVAS 2.1 Introduction A CPVA is composed of a mass that moves relative to the supporting rotor, supported by some type of hinge, typically of a bifilar type using two rollers, as shown in Figure 1.3. Older absorber designs have the absorber center of mass following a circular path [10], but more modern designs use other paths, generally epicycloids, to better handle large amplitude responses [7, 14, 16]. The response of these systems are well captured by models with a degree of freedom for the rotor, and one degree of freedom for each absorber. The governing equations are highly nonlinear and have complicated expressions (involving polynomial and square root terms) resulting from the CPVA path kinematics as it moves with respect to the host rotor. Time integration approaches are easy to apply to this system, yet conducting parameter studies this way requires extensive computational time, since these systems are lightly damped, leading to long transients. Moreover, it is not possible to gather information about unstable responses or possible multiple solutions under given operating conditions using time simulations. Perturbation methods, on the other hand, are quite powerful in terms of obtaining insight about the response characteristics with respect to system parameters, as well as stability characteristics, and these can be applied in some practical parameter ranges. However, they require careful application of parameter scaling and their error bounds are linked to the values of physical parameters that are assumed to be small. In addition, in some practical systems the parameter values do not satisfy the required conditions for using perturbation methods. Also, for systems with absorbers tuned to different orders, this approach can be quite cumbersome [18, 19]. Thus, there is a need for more general and 10 powerful numerical tools for investigating these systems. Here we will elaborate the development of a new numerical solution approach with the intent to circumvent these shortcomings. The approach utilizes a variation of the harmonic balance method (HBM), one of the many standard procedures used in the analysis of nonlinear oscillatory systems. Before going into the details of the particular method discussed here, let us summarize the basic idea of the HBM. The procedure starts by assuming a Fourier series solution truncated to H number of harmonics (plus a DC term, if needed) for the coordinate(s) of the oscillatory system. This series represents a periodic solution corresponding to a steady-state response, without regard to its stability. Thus, the first and second time derivatives of this assumed solution are also readily available up to H harmonics. By plugging these solution forms into the equations of motion and collecting the coefficients of all the harmonics that appear in the resulting expressions, one can obtain algebraic relations between the system parameters and the amplitudes and phases of the response up to the H th harmonic. If this system is under periodic forcing, one would express the truncated Fourier series such that all of the forcing frequencies has a matching harmonic term in the assumed solution. Then the resulting algebraic amplitude and phase relations will depend on the forcing amplitude. If the system under investigation is a free oscillation problem, the base frequency of the Fourier series (the lowest non-DC harmonic) will be a free variable and the resulting algebraic equations will result in an implicit relation between response amplitudes, phases, and frequency. An example to such problem is the determination of response frequency of a simple pendulum as a function of its swing angle [23, 24]. The harmonic balance method presents certain advantages over perturbation methods in the analysis of nonlinear oscillatory systems. First, it does not impose any limitation on the dominance of the nonlinearities in the equations, for example, system parameters linked with the nonlinearities do not have to be small. Secondly, the solutions obtained through HBM are already in the frequency domain, and thus nonlinear interactions are easily investigated. Although the method is quite straightforward, there are several shortcomings of the 11 HBM, in addition to those associated with other purely numerical tools. First, harmonics higher than H appear in the formulation once the assumed solutions are substituted into the nonlinear terms of the equations of motion. These higher harmonics, which carry information on the dynamics of the system, cannot be used in the solution of the response as they are higher than the truncation level H of the assumed solution. Thus, some information about the higher order harmonics is discarded. For example, a cubic nonlinearity would generate harmonics up to 3H. One can generally select H large enough for sufficient accuracy in virtually all cases for the present system, although motions near the cusp point can cause convergence difficulties. The second disadvantage of this method becomes more dominant in the case where the number of harmonics in the assumed solution is increased to alleviate the loss of accuracy discussed above. That is, as H is increased, the complexity of the algebraic expressions generated by the method rapidly increases to the extend that even numerical solution methods would not be able to solve them efficiently. The third problem with this method is that if the nonlinearities in the system involve any terms with non-integer powers, such as square roots, generally it is not possible to reduce the resulting expression in a summation of harmonics through trigonometric identities. The method we use in this study, proposed in [21], eliminates the above mentioned problems related with HBM, to the extend that its usage becomes very suitable to numerically simulate steady state response characteristics of systems involving CPVAs. The main requirement of this method is that the equations of motion of the system being analyzed can be put into a standard form through a series of auxiliary variable definitions. This pushes the information that is normally lost in the HBM operation down to the new algebraic equations that are defined along with the new auxiliary variables. The equations of motion need to satisfy the following form after this recasting process, ˙ = c + l(Z) + q(Z, Z) m(Z) Z ∈ RNeq (2.1) where Z is the vector of Neq unknowns, c is a constant vector, m and l are linear operators, 12 and q is a quadratic operator. As this form has at most quadratic terms, higher harmonics generated by this multiplication are limited to 2H. Moreover, this quadratic form allows for automated generation of the resulting algebraic equations, as well as their gradients, with respect to Z for any truncation level H, whereas trying to apply the HBM directly to the original equations of motion with a high truncation level would result in a set of algebraic equations difficult to manage. Application of this form of the HBM involves expressing Z in a Fourier series as H Z(t) = Z0 + H Zc,k cos(kωt) + k=1 Zs,k sin(kωt) (2.2) k=1 and substituting this form into Equation (2.1). Solving the resulting algebraic system of equations for the Fourier coefficients of each harmonic provides an approximate form for a periodic solution to the original system. In order to systematically obtain the algebraic equations, the coefficients are collected in a single vector U = Z0 , Zc,1 , Zs,1 , Zc,2 , Zs,2 , . . . , Zc,H , Zs,H . (2.3) The algebraic system of coefficients then takes the form ωM (U ) = C + L(U ) + Q(U, U ) (2.4) where the operators M (.), C, L(.) and Q(U, U ) depend only on the original operators m(.), c, l(.), q(., .) and the number of harmonics assumed in Equation (2.2), H. The details of the relations between these sets of operators are derived in [21]. This means, as long as the original system can be put into form Equation (2.1), the system of equations resulting from the application of the HBM can be automatically generated. Numerical solutions to the resulting system of quadratic algebraic equations can be obtained using one of many available algorithms. The algorithm used in [21] is a pseudoarclength continuation method, codes of which are freely available to public under the package name called MANLAB [25], where one of the system parameters is chosen as a free 13 variable and periodic solutions are obtained by solving the Fourier coefficients as a function of this parameter. Moreover, the stability information of these periodic solutions can be obtained through the methodology used in [26], which first assembles the Hill’s matrix truncated at H th harmonics at any given solution point. This matrix contains individual harmonics of the Jacobian of the original equations of motion as block sub-matrices and will be different at any solution point, the stability of which being investigated. The algorithm then uses the eigenvalues of the resulting Hill’s matrix to calculate the Floquet exponents of the corresponding periodic solution in order to assess its stability. Expressing the Jacobian of the equations of motion in terms of the solutions found using HBM can be quite difficult. In this work this process is facilitated by the following approach: we reconstruct the calculated periodic solutions in the time domain, substitute these into the Jacobian expression, and use a Fast Fourier Transform to obtain the individual harmonics of the Jacobian. In this chapter, we describe the application of this methodology to a quite generic system model for a rotor fitted with multiple CPVAs. The model consists of a rigid rotor under single harmonic order excitation and N point mass CPVAs fitted to this rotor and allowed to move along paths described by Section 1.1. Several other configurations, such as CPVAs under the influence of gravity, CPVAs with rocking motion (studied in Chapter (4)) were also successfully analyzed using this method. Although similar in essence, this simplest form is chosen for presentation here in order to demonstrate the essential steps in the development and implementation of this approach. 2.2 Dynamical Model The system illustrated in Figure 2.1, incorporates point mass CPVAs prescribed to move along indicated paths with respect to the rotor. The excitation torque applied to the rotor has order-domain characteristics, namely the frequency of the excitation is a function of 14 θ R(S) ρ0 S Figure 2.1: Schematic of a rotor fitted with CPVAs. rotor angular position, θ(t). This system is modeled using Lagrange’s equations for a rigid rotor with N CPVAs, resulting in a system with N + 1 degrees of freedom. The equations of motion of the system are given by J θ¨ = T (θ) + T0 − c0 θ˙ N − mj Rj2 (Sj )θ¨ + Gj (Sj )S¨j j=1 2 dGj (Sj ) 2 dRj (Sj ) ˙ ˙ Sj + S˙ j θ + dSj dSj dRj2 (Sj ) 1 2 ¨ ˙ ¨ mj Sj + mj Gj (Sj )θ + mj θ = −ca,j S˙j 2 dSj (2.5) (2.6) where j = 1 : N where Sj is the (arc length) displacement of the j th absorber COM relative to the vertex on its path, mj is the mass of absorber j, Rj (Sj ) and Gj (Sj ) are functions related to the path of the absorber as defined in Section 1.1, J is the moment of inertia of the rotor, T (θ) is the applied torque acting on the rotor, ca,j represents the equivalent viscous damping coefficient for absorber j as it moves along its path, T0 is the constant torque applied to the rotor that maintains a mean operating rotational speed of Ω by counteracting viscous rotor damping that is modeled by coefficient c0 . 15 The forcing term T (θ) in Equation (2.5) has the form of engine order excitation, so that rather than explicitly depending on time it is a function of rotor position. It is reasonable to assume that the rotor angle θ(t) can be used as the independent variable in the forcing term, as T (t) = T sin(nθ) (2.7) where n is the engine excitation order. For example, in four-stroke piston engines, n is equal to one half the number of cylinders. 2.2.1 Non-dimensionalization Table 2.1: Definition of dimensionless variables. Variables Definitions Variables Definitions αj mj Rj2 (0)/J gj (sj ) Gj (Sj )/Rj (0) µa,j ca,j/mj Ω rj2 (sj ) Rj2 (Sj )/R2 (0) j µ0 c0/JΩ Γ T/JΩ2 sj Sj/Rj (0) Γ0 T0/JΩ2 The set of equations are non-dimensionalized using the definitions listed in Table 2.1. As time does not appear explicitly in the equations, it is convenient to replace time as the independent variable with the rotor position θ. This results in the following definition for the dimensionless speed of the rotor, ν(θ) = θ˙ . Ω (2.8) These operations and definitions, along with use of the chain rule, yield the following dimensionless set of equations that govern the behavior of the rotor speed and the absorber 16 displacements d(rj2 (sj )) N νν = − αj gj (sj ) νν sj + ν 2s j + rj2 (sj )νν j=1 − νsj gj (sj )µa,j νsj + ν + ν 2s j dsj + ν 2 sj2 + Γ sin(nθ) + Γ0 − µ0 ν 2 1 d(rj (sj )) sj + gj (sj ) − ν = −sj µa,j 2 dsj d(gj (sj )) dsj (2.9) (2.10) d() where () = dθ . Note that these equations now have the familiar form of periodic excitation, where nθ plays the role usually taken by ωt. These equations have been the basis of many previous investigations of CPVA systems, including [16, 15, 19, 18, 20, 7, 8]. 2.2.2 Reforming the equations for HBM In order to put the non-dimensional equations of motion given in Equations (2.9 and 2.10) into the form described in Equation (2.1), a large number of new variable definitions are required. The procedure starts by inserting the coordinates and their derivatives of the original equations into the new unknown vector, Z. z1 = ν (2.11) z2 = ν (2.12) z(3j) = sj (2.13) z(3j+1) = sj (2.14) z(3j+2) = sj . (2.15) 17 The next step is to define variables to represent the radial positions of the absorber COM and their first and second derivatives, as z(3N +3j) = rj2 (sj ) z(3N +3j+1) = z(3N +3j+2) = (2.16) d(rj2 (sj )) (2.17) dsj d2 (rj2 (sj )) (2.18) ds2j such that 2 4 + φ1,j z(3j) z(3N +3j) = 1 − n ˜ 2j z(3j) (2.19) 3 z(3N +3j+1) = −2˜ n2j z(3j) + 4φ1,j z(3j) (2.20) 2 . z(3N +3j+2) = −2˜ n2j + 12φ1,j .z(3j) (2.21) In order to account for powers of the absorber position and velocity coordinates that are higher than quadratic, another set of variables are defined as z(6N +2j+1) = s2j 2 = z(3j) (2.22) z(6N +2j+2) = sj2 2 = z(3j+1) . (2.23) The most problematic terms in the equations of motion stem from the functions gj (sj ) and their derivatives d(gj (sj )) , sj as the definition of the former involves the square root of an expression (see Equation (1.2)). To handle this, we define a new set of variables and form a set of quadratic algebraic equations as follows z(8N +2j+1) = gj (sj ) = 1 rj2 (sj ) − 4 d(gj (sj )) 1 z(8N +2j+2) = = dsj 2gj (sj ) d(rj2 (sj )) 2 dsj d(rj2 (sj )) dsj 2 2 1 d (rj (sj )) − 2 ds2j (2.24) (2.25) to obtain 1 2 2 z(8N +2j+1) = z(3N +3j) − 4 z(3N +3j+1) 1 2z(8N +2j+2) z(8N +2j+1) = z(3N +3j) − z(3N +3j+1) z(3N +3j+2) . 2 18 (2.26) (2.27) Up to this point, the above operations yield the following form z(1) = z(2) (2.28) z(3j) = z(3j+1) (2.29) z(3j+1) = z(3j+2) (2.30) N 0 = −Γ0 − Γ sin(nθ) + µ0 z(1) + z(1) z(2) + z(1) αj − z(3j+1) z(8N +2j+1) µa,j j=1 + z(1) z(3j+1) z(3N +3j+1) + z(3j+2) z(8N +2j+1) + z(6N +2j+2) z(8N +2j+2) + z(2) z(3j+1) z(8N +2j+1) + z(3N +3j) (2.31) 1 0 = µa,j z(3j+1) + z(1) z(3j+2) − z(3N +3j+1) + z(2) z(3j+1) + z(8N +2j+1) 2 (2.32) 2 0 = 1 − z(3N +3j) − n ˜ 2j z(6N +2j+1) + φ1,j z(6N +2j+1) (2.33) 0 = −z(3N +3j+1) − 2˜ n2j z(3j) + 4φ1,j z(3j) z(6N +2j+1) (2.34) 0 = −z(3N +3j+2) − 2˜ n2j + 12φ1,j z(6N +2j+1) (2.35) 2 0 = −z(6N +2j+1) + z(3j) (2.36) 2 0 = −z(6N +2j+2) + z(3j+1) (2.37) 1 2 2 0 = z(3N +3j) − z(3N +3j+1) − z(8N +2j+1) 4 1 0 = z(3N +3j+1) − z(3N +3j+1) z(3N +3j+2) − 2z(8N +2j+1) z(8N +2j+2) 2 (2.38) (2.39) The only terms that do not satisfy the required form are left in Equation (2.31), which governs the rotor dynamics. To handle these, we define N z(10N +3) = αj z(3j+1) z(3N +3j+1) + z(3j+2) z(8N +2j+1) + z(6N +2j+2) z(8N +2j+2) j=1 (2.40) N z(10N +4) = αj z(3N +3j) + z(3j+1) z(8N +2j+1) j=1 19 (2.41) N z(10N +5) = − αj µa,j z(3j+1) z(8N +2j+1) + z(1) z(10N +3) + z(2) z(10N +4) . (2.42) j=1 It is also convenient to have the dimensionless rotor acceleration (νν ) as one of the components of the solution, since it is an important output of the response, needed to assess the vibration reduction of the rotor. So, we define it as z(10N +6) = z(1) z(2) . (2.43) This completes the recasting operation and yields the desired set differential-algebraic set of equations in their nearly final form. What remains to be done is to deal with the two excitation terms, Γ0 and Γ sin(nθ). Γ0 is the dimensionless mean torque used to maintain a steady mean speed (ν ≈ 1), and thus we take Γ0 = µ0 , which balances the DC components to leading order. Then, after the sine and cosine harmonic amplitudes of the unknown vector Z have been assembled, the harmonic excitation torque amplitude Γ is added to the algebraic equation in the set Equation (2.4) corresponding the harmonic sin(nθ). This yields the final form of the equations to be solved: z(1) = z(2) (2.44) z(3j) = z(3j+1) (2.45) z(3j+1) = z(3j+2) (2.46) 0 = −Γ0 − Γ sin(nθ) + µ0 z(1) + z(1) z(2) + z(1) z(10N +5) (2.47) 1 0 = µa,j z(3j+1) + z(1) z(3j+2) − z(3N +3j+1) + z(2) z(3j+1) + z(8N +2j+1) 2 (2.48) 20 2 0 = 1 − z(3N +3j) − n ˜ 2j z(6N +2j+1) + φ1,j z(6N +2j+1) (2.49) 0 = −z(3N +3j+1) − 2˜ n2j z(3j) + 4φ1,j z(3j) z(6N +2j+1) (2.50) 0 = −z(3N +3j+2) − 2˜ n2j + 12φ1,j z(6N +2j+1) (2.51) 2 0 = −z(6N +2j+1) + z(3j) (2.52) 2 0 = −z(6N +2j+2) + z(3j+1) (2.53) 1 2 2 0 = z(3N +3j) − z(3N +3j+1) − z(8N +2j+1) 4 1 0 = z(3N +3j+1) − z(3N +3j+1) z(3N +3j+2) − 2z(8N +2j+1) z(8N +2j+2) 2 (2.54) (2.55) N 0 = −z(10N +3) + αj z(3j+1) z(3N +3j+1) + z(3j+2) z(8N +2j+1) j=1 + z(6N +2j+2) z(8N +2j+2) (2.56) N 0 = −z(10N +4) + αj z(3N +3j) + z(3j+1) z(8N +2j+1) (2.57) αj µa,j z(3j+1) z(8N +2j+1) + z(1) z(10N +3) + z(2) z(10N +4) (2.58) j=1 N 0 = −z(10N +5) − j=1 0 = −z(10N +6) + z(1) z(2) (2.59) As noted above, standard numerical tools are used to solve these equations as a single parameter is varied. For CPVA systems, the most natural parameter to vary is the fluctuating torque amplitude Γ, since in a physical system all other parameters will be fixed and the system operation is evaluated by considering the operating range of Γ. 2.3 Case Studies In order to test the capabilities of this solution approach, certain configurations of the system that are known to produce rich dynamical behavior are investigated. Perturbation methods such as the method of averaging can accurately predict the dynamical behavior of these systems, provided assumptions hold about small values of certain non-dimensional parameters 21 [27, 28, 29, 15, 20] . Specifically, it is common to choose the ratio of the total rotational inertia of the absorbers to the rotor inertia as the main small parameter, since it is typically less than 0.15 in practice. Other parameters, such as the absorber damping and its detuning away from the engine order, and the excitation torque amplitude, must also be small when scaled appropriately. However, in systems with an inertia ratio that is not small, one generally relies on time integration of the equations of motion, which, as mentioned before, are time consuming and cumbersome for carrying out parameter studies. Also, perturbation methods assume a single harmonic response of the absorbers, which is a good approximation in many cases, but is known to break down, for example, when the absorbers move at amplitudes near the cusp of epicycloidal and cycloidal paths. Here, we show case studies for systems where the parameter values do not permit the use of perturbation methods, and also generate solutions that capture several harmonics. Furthermore, the harmonic balance approach captures unstable responses that cannot be found by time simulations. 2.3.1 Response of Circular Path CPVA with High Inertia Ratio The first case considered is a rotor equipped with a single circular path CPVA whose inertia is not small compared with that of the rotor. It is well known [11, 12] that circular path CPVAs tuned close to the excitation order will undergo a jump instability at a critical excitation torque amplitude. These jumps do not occur in CPVAs with paths such as cycloids or the tautochronic epicycloids, unless they are undertuned [7]. This renders the use of circular path CPVAs impractical at close tunings to the excitation order, as this jump instability results in a phase shift that causes the CPVA to amplify the rotor torsional vibrations. Therefore, such absorbers are generally overtuned, which reduces their effectiveness but extends their operating range. In the present model, we consider a single harmonic engine excitation of order n = 2 and a rotor/CPVA system for which the inertia ratio is α = 0.9, which is far too large for the application of perturbation methods. The steady-state response characteristics of this system as a function of the excitation torque amplitude are given 22 in Figure 2.2, for both the dimensionless rotor acceleration and the absorber displacement response. Amplitudes and phases for several harmonics are shown. Of course, the total response is a combination of these harmonics, and the response waveform depends on the relative phases of the harmonics. The response waveforms of several solution points are presented in Figures 2.4 - 2.7. The overall response as a function of the torque amplitude Γ has the classic form for a softening nonlinear system, that is, as the torque increases the amplitude increases until it hits a turning point (a saddle-node bifurcation), and the there exists another, larger amplitude saddle-node bifurcation, such that there are three branches of response over a range of torques, the upper and lower of which are generally dynamically stable and the middle of which is unstable. (For systems with multiple CPVAs the picture is more complicated, as shown in the next example.) A stability analysis of these periodic solutions is conducted using the Floquet analysis described above. The indicator of stability of a given solution is the amplitude of the maximum Floquet multiplier of that solution obtained through the algorithm. These maximum Floquet multipliers are plotted against Γ in Figure 2.3. On the lower branch, the order n absorber phase is approximately π, and thus it counteracts the applied torque, resulting in vibration reduction. On the upper response branch the CPVA acts like a vibration amplifier, so this response must be avoided in practice. As can be seen from the harmonic amplitudes, on the lower response branch the absorber is dominated by the n = 2 harmonic, although other harmonics contribute significantly to the response on the upper branch. On the other hand, the harmonic content of the rotor response is dominated by higher order harmonics in all cases, since the n = 2 harmonic is very small, as it should be when the absorber is operating to reduce the order n rotor response. Overall, the higher harmonic amplitudes do not sequentially diminish, for instance, the order n = 6 harmonic amplitude of the rotor response is larger than the order n = 4 amplitude for certain torque ranges. 23 0.35 0.35 Order n Locked Abs. Order n Order 2 n Order 3 n 0.25 0.2 0.15 0.1 0.05 0.25 0.2 0.15 0.1 0.05 0 0 0 0.1 0.2 0.3 0.4 Excitation Torque (Γ) 0.5 0 (a) Rotor Acceleration Amplitudes 0.1 0.2 0.3 0.4 Excitation Torque (Γ) 0.5 (b) Absorber Displacement Amplitudes 200 200 150 150 Absorber Response Phase (deg) Rotor Acceleration Phase (deg) CPVA #1, Order 1 n CPVA #1, Order 2 n CPVA #1, Order 3 n 0.3 Absorber Response Mag. (|s|) Rotor Acceleration Mag. (|ν ν'|) 0.3 100 50 0 -50 -100 -150 100 50 0 -50 -100 -150 -200 -200 0 0.1 0.2 0.3 0.4 Excitation Torque (Γ) 0.5 0 (c) Rotor Acceleration Phases 0.1 0.2 0.3 0.4 Excitation Torque (Γ) 0.5 (d) Absorber Displacement Phases Figure 2.2: Harmonic amplitudes of the response of a circular path CPVA system, where N = 1, n ˜ 1 = 2.01, α1 = 0.9, µa,1 = 0.005, H = 48, n = 2. Only the n, 2n, and 3n harmonics, out of 48, are shown. The black curve in (a) indicates the reference rotor response at order n when the absorber is locked. Dashed parts of the curves indicate unstable periodic solutions. 24 Max Floquet Multiplier 10 2 10 1 10 0 10 -1 0 0.1 0.2 0.3 0.4 Excitation Torque (Γ) 0.5 Figure 2.3: Maximum Floquet multiplier magnitudes of the system in Figure 2.2 at corresponding periodic solutions. Values above one indicate an unstable solution. 25 Excitation Torque (Γ) = 0.45 0.3 Rotor Acceleration (ν ν') 0.2 0.1 0 -0.1 -0.2 -0.3 -0.4 0 1 2 3 4 5 6 5 6 Rotor Revolutions (a) Rotor Acceleration Excitation Torque (Γ) = 0.45 0.2 Absorber Displacement (s) 0.15 0.1 0.05 0 -0.05 -0.1 -0.15 -0.2 0 1 2 3 4 Rotor Revolutions (b) Absorber Displacement Figure 2.4: Rotor angle (acting similar to time) domain reconstruction of a periodic solution of the system in Figure 2.2. This solution belongs the lower stable branch at Γ = 0.45. The blue curve is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. (a) Rotor Acceleration, (b) Absorber Displacement 26 Excitation Torque (Γ) = 0.37 4 Rotor Acceleration (ν ν') 3 2 1 0 -1 -2 0 1 2 3 4 5 6 5 6 Rotor Revolutions (a) Rotor Acceleration Excitation Torque (Γ) = 0.37 0.6 0.5 Absorber Displacement (s) 0.4 0.3 0.2 0.1 0 -0.1 -0.2 -0.3 -0.4 0 1 2 3 4 Rotor Revolutions (b) Absorber Displacement Figure 2.5: Rotor angle domain reconstruction of a periodic solution of the system in Figure 2.2. This solution belongs the middle unstable branch at Γ = 0.37. The blue curve is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. Note that even numerical solution and integration tolerances are sufficient to cause an ultimate deviation from the periodic solution. In this case, the absorber response grew to a point that the numerical integration crashed. (a) Rotor Acceleration, (b) Absorber Displacement 27 Excitation Torque (Γ) = 0.20 0.8 Rotor Acceleration (ν ν') 0.6 0.4 0.2 0 -0.2 -0.4 -0.6 -0.8 0 1 2 3 4 5 6 5 6 Rotor Revolutions (a) Rotor Acceleration Excitation Torque (Γ) = 0.20 0.3 Absorber Displacement (s) 0.2 0.1 0 -0.1 -0.2 -0.3 0 1 2 3 4 Rotor Revolutions (b) Absorber Displacement Figure 2.6: Rotor angle domain reconstruction of a periodic solution of the system in Figure 2.2. This solution belongs the middle unstable branch at Γ = 0.20. The blue curve is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. In this case, the numerically integrated solution eventually jumps down to the lower stable branch. (a) Rotor Acceleration, (b) Absorber Displacement 28 Excitation Torque (Γ) = 0.15 0.8 Rotor Acceleration (ν ν') 0.6 0.4 0.2 0 -0.2 -0.4 -0.6 -0.8 0 1 2 3 4 5 6 5 6 Rotor Revolutions (a) Rotor Acceleration Excitation Torque (Γ) = 0.15 0.4 Absorber Displacement (s) 0.3 0.2 0.1 0 -0.1 -0.2 -0.3 -0.4 0 1 2 3 4 Rotor Revolutions (b) Absorber Displacement Figure 2.7: Rotor angle domain reconstruction of a periodic solution of the system in Figure 2.2. This solution belongs the upper stable branch at Γ = 0.15. The blue curve is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. (a) Rotor Acceleration, (b) Absorber Displacement 29 2.3.2 Bifurcation to Nonunison Response of Two Identical Tautochronic CPVAs In general, balancing and space restrictions in applications force designers to distribute the absorber inertia into several absorbers placed around and along the rotor. The response characteristics of such groups of CPVAs are of importance, since the maximum torque range can be achieved when the CPVAs move in a synchronous manner, that is, all absorbers respond with equal amplitude and phase. In [7], it was shown that sets of multiple identical CPVAs with tautochronic paths undergoing a synchronous response can undergo an instability to a non-synchronous response at a certain torque amplitude that depends on the system parameters. (These tautochronic paths render the absorber essentially linear out to large amplitudes, such that the dominant system nonlinearity comes from the interactions between the rotor and the absorbers, rather than from the absorber path, as in the previous case.) To consider systems of this type, we consider the response of a rotor subjected to order n = 1.5 torque, fitted with two identical absorbers with tautochronic paths, and a small inertia ratio so that the perturbation results of [7] can be verified. The purpose of this case study is to check the ability of the HBM to capture bifurcations to non-synchronous responses, and to track the resulting non-synchronous responses. The solution branches for the amplitudes of several harmonics obtained by the HBM are presented in Figure 2.8, for both the dimensionless rotor acceleration and the normalized absorber displacement. At low torques the two absorbers respond in a synchronous manner and the rotor responds to this motion as if the two absorbers were a single absorber with the same total mass, as can be identified in the plots. As the torque is increased, there is a splitting of the two absorber responses, where the synchronous response becomes unstable and separate branches of absorber response emerge, corresponding to distinct responses of each absorber. Note that this response does not result in a significant increase in the rotor response, but it causes one absorber to reach its cusp (maximum possible response amplitude) at a lower torque than the corresponding synchronous response. This leaves the available amplitude range of the 30 second absorber unused, thereby limiting the feasible operating torque range of the entire system. It is interesting to note that the bifurcation predicted by the perturbation analysis is sharp, whereas the HBM shows a more initially gradual separation of the branches followed by a more dramatic separation. This is due to convergence issues in the HBM solvers, but in fact represents more closely what one would expect to see in practice, where slight differences between the absorbers will result in a similar effect. 31 Rotor Acceleration Mag. (|ν ν'|) 0.15 Order n Locked Abs. Order n 0.1 0.05 0 0.08 0.1 0.12 0.14 Excitation Torque (Γ) 0.16 0.18 0.2 (a) Rotor Acceleration Amplitudes Absorber Response Mag. (|s/smax |) 1 0.9 0.8 0.7 0.6 0.5 CPVA #1, Order n CPVA #2, Order n 0.4 0.08 0.1 0.12 0.14 Excitation Torque (Γ) 0.16 0.18 0.2 (b) Absorber Displacement Amplitudes Figure 2.8: Harmonic amplitudes of the response of the rotor with two identical tautochronic CPVAs, with N = 2, n ˜ 1,2 = 1.50, α1,2 = 0.1, µa,1,2 = 0.005, H = 16, n = 1.5 and the order of the first assumed harmonic is 0.5. (a) Rotor Acceleration, (b) Absorber Displacement 32 Excitation Torque (Γ) = 0.12 0.06 Rotor Acceleration (ν ν') 0.04 0.02 0 -0.02 -0.04 -0.06 -0.08 0 0.5 1 1.5 2 2.5 3 Rotor Revolutions (a) Rotor Acceleration Excitation Torque (Γ) = 0.12 0.3 Absorber Displacement (s) 0.2 0.1 0 CPVA #1 CPVA #2 -0.1 -0.2 -0.3 0 0.5 1 1.5 2 2.5 3 Rotor Revolutions (b) Absorber Displacement Figure 2.9: Rotor angle domain reconstruction of a periodic solution of the system in Figure 2.8. The excitation torque Γ = 0.12 in this case correspond to a solution where the two CPVAs move in unison. The blue curve in (a) is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. Numerical integration results of the CPVA signals are indistinguishable from the HBM results and are not shown for the sake of clarity in the figure. (a) Rotor Acceleration, (b) Absorber Displacement 33 Excitation Torque (Γ) = 0.14 0.1 0.08 Rotor Acceleration (ν ν') 0.06 0.04 0.02 0 -0.02 -0.04 -0.06 -0.08 0 0.5 1 1.5 2 2.5 3 2.5 3 Rotor Revolutions (a) Rotor Acceleration Excitation Torque (Γ) = 0.14 0.4 Absorber Displacement (s) 0.3 0.2 0.1 CPVA #1 CPVA #2 0 -0.1 -0.2 -0.3 -0.4 0 0.5 1 1.5 2 Rotor Revolutions (b) Absorber Displacement Figure 2.10: Rotor angle domain reconstruction of a periodic solution of the system in Figure 2.8. This solution is where the two CPVAs are non-synchronous at Γ = 0.14. The blue curve in (a) is the signal generated by evaluating the Fourier coefficients of the HBM solution, and the red dashed curve is the result of numerical integration started at the same initial conditions as the beginning of the blue signal. As above, numerical integration results of the CPVA signals are not shown. (a) Rotor Acceleration, (b) Absorber Displacement 34 2.4 Conclusion In this chapter we have presented the application of a harmonic balance method for determining the steady-state response of a rotor systems equipped with CPVAs. The nonlinear equations of motion for these systems are generally analyzed either with time integration or using perturbation methods, both of which have limitations in terms of speed or applicability. The current method is adopted from [21], in which the nonlinearities are put into a quadratic polynomial form. The solution is expressed as a Fourier series for application of the HBM with a relatively high number of assumed harmonics. Essentially the method differs from the classical harmonic balance method in that it has a pre-processor for the problem formulation. It operates on the system equations of motion and rewrites them as a larger, yet equivalent, set of equations using auxiliary variables. In this way, the application of the HBM, namely assuming a periodic solution with a given number of harmonics and constructing an algebraic system of equations that govern the amplitudes of the harmonics of the dependent variables, can be automated irrespective of the structure of the original equations of motion. Estimation of the stability of the periodic steady state solutions obtained in this way was also made possible using a related algorithm outlined in [26]. The requirement for this secondary step in the solution process is that each component of the Jacobian of the original equations of motion must be written in terms of its harmonic coefficients at each periodic solution, where stability information is sought. The form of the Jacobian of systems with CPVAs is much more complex than the equations of motion themselves. Therefore, here we evaluated it numerically by reconstructing the periodic solutions of the states of the equations of motion using the coefficients obtained by HBM and extracted the Jacobian harmonic coefficients numerically using a Fast Fourier Transform. This solution framework was put to test for two different system configurations that are known to show significant qualitative changes in their response characteristic as the 35 excitation amplitude is varied. The results showed that this approach can effectively be used in the analysis of these CPVA systems. The complexity of the CPVA systems investigated here was kept at a level that makes the demonstration of the method as clear as possible. However, more involved configurations have been successfully analyzed with this algorithm. This tool will allow for the rapid study of absorber systems with parameter values that render analytical tools inapplicable and/or direct numerical simulations too slow. Such a tool will be very useful for parameter studies and design purposes. 36 CHAPTER 3 LINEAR AND NONLINEAR DYNAMICS OF CPVA SYSTEMS ON FLEXIBLE ROTOR SYSTEMS 3.1 Introduction Flexible rotating shafts, such as crankshafts of internal combustion engines, are subject to engine order excitation, typically composed of several harmonics. Due to the fact that an engine order excitation frequency is proportional to the engine speed, certain harmonics of this excitation can coincide with the resonant frequencies of the flexible shaft as the rotor speed varies, resulting in problematic resonant vibrations which can lead to failure Figure 3.1. Torsional versions of tuned mass absorbers (TMA) are commonly used to rectify the adverse effect of these possible resonance conditions [30, 31, 32]. A TMA is tuned to the problematic resonance frequency and the adapted system has two resonances surrounding the original untreated resonance. The damping of the TMA limits the new resonance amplitudes, resulting in reduced vibration levels across the critical speed. Such TMAs are used to address both flexural and torsional shaft modes. In contrast, CPVAs are torsional vibration reduction elements widely used in light aircraft engines, and being considered for modern automotive powertrains as well, where they allow engines to run in efficient low-speed operation with minimal vibration [9, 33]. Unlike a TMA, CPVAs do not incorporate elastic elements, and their effective stiffness stems from the centrifugal field. Therefore, they can be tuned to address a certain engine order over all engine speeds, whereas TMAs are tuned to a certain frequency and are effective over a limited speed range. Also, CPVAs are generally used to address torsional vibrations, but can be adapted for more general use [34, 35]. In much of the previous work on CPVAs a rigid rotor assumption has been used, and the primary intention of the CPVA is to reduce torsional vibration of a rotating shaft that 37 Figure 3.1: A crankshaft with torsional failure [3] arises from a dominant harmonic of engine order excitation [29, 16]. Such analyses are valid as long as the product of engine order and operating speed remains well below the natural frequencies of any rotor torsional modes. In this chapter, we consider the response of a model for a flexible rotor fitted with CPVAs. Linear modal analysis that incorporates rotor torsional vibration shows that CPVAs can also be used to modify the structure of resonant frequencies of the coupled rotor-absorber system, due to eigenfrequency veering behavior. In fact, the analysis indicates that certain resonances can be completely eliminated for all possible engine speeds [17, 22, 36]. The key plots and ideas that show this behavior are given in the subsequent sections. While this linear vibration analysis offers significant insight about the system dynamics, in practice nonlinear effects from the CPVAs must be taken into account, since these devices are generally designed to operate over large amplitudes. We consider the simplest model that includes the features of interest, namely a torsional vibration mode, a CPVA, and engine order excitation. Thus, the rotor is modeled using two lumped rotational inertial elements connected with a linear torsional spring, so that it has a speed-independent torsional natural frequency before the addition of the CPVA. In the model, a CPVA is added to one rotor element while an engine-order excitation torque is applied to the other. The CPVA is modeled as a point mass that travels along a prescribed 38 path relative to its host rotor. The equations of motion for this model are obtained using Lagrange’s method, and the equations are non-dimensionalized and scaled with the assumptions that the effective inertia of the CPVA is much smaller than that of the rotor system. It is also assumed that the non-dimensional versions of the amplitude of the excitation torque, the absorber damping, the relative deflection between the rotors, and the deviation of the CPVA path from a tautochronic epicyloid curve are also small [7, 14, 19]. These assumptions are consistent with practical absorber implementations and make the system amenable to analytical treatment. These assumptions lead to a set of two weakly coupled, weakly nonlinear, oscillator equations driven near resonance, which govern the absorber dynamics and the rotor torsional deflection. Near the speed range where the frequency veering occurs these equations have a one-to-one internal resonance, and this is precisely the parameter region of primary interest in applications. The method of averaging is applied to these equations to investigate the system response as function of the system and excitation parameters. The model shows that the favorable features of the linearized system model are retained when the absorber amplitudes become large, so long as the absorber path is taken to be close to a tautochronic epicycloid. This is crucial for designs in which the absorber mass and/or placement is constrained, or engine speeds are low, resulting in large amplitude motions of the absorbers. The remainder of this chapter is organized with the following sections: the dynamical model, which outlines the assumptions and derivation of the equations of motion, and provides definitions for dimensionless parameters and the path followed by the absorber mass; the system response, which describes the scaling of parameters and the application of the method of averaging; parameter studies, which provides a sampling of responses obtained from the averaged equations; and conclusions, which summarizes the results and their significance. 39 3.2 Dynamical Model S θ1 T θ2 m R(S) k1 J J 1 2 Figure 3.2: Illustration of the model for a single rotor torsional mode and a single CPVA. The system shown in Figure 3.2 and described above is modeled using Lagrange’s equations, and the equations of motion are given by J1 θ¨1 + k1 (θ1 − θ2 ) = T (t) d R2 (S) 1 m S¨ + θ¨2 G (S) − mθ˙22 = −ca S˙ 2 dS d (G (S)) ˙ 2 J2 + m(R2 (S)) θ¨2 + mG (S) S¨ + m S dS d R2 (S) ˙ (S) + k1 (θ2 − θ1 ) = ca SG +mθ˙2 S˙ dS (3.1) (3.2) (3.3) where S is the displacement of the absorber mass m relative to its vertex, θ1,2 represent the angles of the rotor elements, R(S) and G(S) are functions related to the path of the absorber as described in Section 1.1, J1,2 are the moments of inertia of the rotor elements, k1 is the torsional stiffness of the inter-rotor torsional spring, T (t) is the applied torque acting on the first rotor element (described in more detail below), and ca represents the equivalent viscous damping coefficient for the absorber as it moves along its path. Note that the bearing resistance for the rotors is not included in this model, nor is the mean (DC) component of the applied torque, but generally these balance to set the mean rotor speed Ω. 40 Of primary interest in applications are the rigid body torsional vibration and the relative deflection between the rotors. To this end we switch to coordinates expressed in terms of the sum and difference of the absolute angle coordinates of the two rotors, defined as θ (t) + θ2 (t) θ (t) − θ1 (t) z(t) = 1 , w(t) = 2 2 2 (3.4) which results in equations of motion of the form J1 (¨ z − w) ¨ − 2k1 w = T (t) d R2 (S) 1 m S¨ + G(s) (¨ z + w) ¨ − m (z˙ + w) ˙ 2 = −ca S˙ 2 dS (3.5) (3.6) J2 + m(R2 (s)) (¨ z + w) ¨ + k1 z = −m G(S)S¨ + d (G (S)) ˙ 2 d R2 (S) ˙ S + S (z˙ + w) ˙ dS dS ˙ + ca G(s)S. (3.7) The forcing term T (t) in Equation (3.1) has the form of engine order excitation, so that rather than explicitly depending on time it is a function of rotor position. It is reasonable to assume that the mean rotor angle z(t) can be used as the independent variable in the forcing term, since relative deflections of the rotor will be small, so that the dominant harmonic of the torque can be expressed as T (z) = T sin(nz). 3.2.1 (3.8) Non-dimensionalization The set of equations are non-dimensionalized using the definitions listed in Table 3.1. As none of the terms explicitly depend on time in the equations, it is convenient to replace time as the independent variable with the mean rotor position z. This results in the following definition for the dimensionless mean speed of the average rotor position: ν(z) = 41 z˙ . σ (3.9) Table 3.1: Definition of the non-dimensional variables. Variables Definitions (J1 +J2 )k1 J1 J2 S R(0) G(S) R(0) R2 (S) R2 (0) ω0 s g(s) r2 (s) T J1 Ω 2 mR2 (0) J2 Ω ω0 Γ α σ These operations and definitions, along with the simplifying assumption J1 = J2 , yield the following dimensionless set of equations that govern the behavior of the rotor rigid body mode, the relative rotor deflection, and the absorber: 1 νν = − αν ν g(s)s + r2 (s) ν (1 + w ) + νw 2 + ν g(s)s + s s dg(s) d(r2 (s)) + (1 + w ) ds ds + Γ sin(nz) 2 2 1 2 d(r (s)) νs + νg(s)w = ν 1 + w + ν s + g(s)(1 + w ) − µa s 2 ds 1 ν 2 σ 2 w + w = −σ 2 ν ν w + α r2 (s) ν + ν w + νw 2 dg(s) + µa s g(s) + ν s g(s) + νs g(s) + νs s ds 2 d(r (s)) + νs 1 + w , ds d() where () = dz . This dynamic model is the basis for the analysis that follows. 42 (3.10) (3.11) (3.12) 3.3 System Response The behavior of the natural frequencies of the linearized system is the primary motivating factor for extending the analysis to include nonlinear behavior. The inclusion of a CPVA into a flexible rotor system results in an eigenfrequency veering behavior as the mean rotor speed is varied. This results in a gap between the two nonzero eigenfrequencies of this three degree-of-freedom system (one of the eigenfrequencies is always zero due to the rigid body motion). The width of the gap depends on the inertia ratio parameter α, which acts as the coupling parameter in this system. To illustrate this, Equations (3.10 - 3.12) are linearized and the eigenfrequencies are plotted as a function of the dimensionless rotor speed σ for α → 0 and for α = 0.005 and α = 0.01 in Figure 3.3 . Note that the independent variable in Equations (3.10 - 3.12) is the mean rotor position z, and therefore the eigenfrequencies obtained are in the (nondimensional) order domain, rather than the frequency domain, as can be seen in Figure 3.3 . Moreover, a given excitation order corresponds to a horizontal line in these plots. If the system parameters are tuned such that a problematic order excitation line passes through the veering gap, both resonances are avoided [36]. This is an extremely intriguing possibility, but if this is to be a feasible tuning strategy in practice, the nonlinear system response must be considered, to which we turn next. 3.3.1 Scaling and Averaging A scaling scheme is applied to the equations of motion so that the contributions of the terms that are small in practice fall into the right order for the purpose of applying a perturbation analysis. To this end, the following scaling definitions are used in Equations (3.10 - 3.12), where 1 is the central small parameter, ˆ µa = µ α = 2 , w = w, ˆ Γ = 2 Γ, ˆa , φ = φˆ (3.13) which capture small values for the dimensionless parameters associated with the inertia ratio, the absorber damping, the relative rotor deflection, the applied torque, and the deviation 43 2.6 Ord e r 2.4 2.2 2 1.8 1.6 0.4 0.45 0.5 σ 0.55 0.6 Figure 3.3: Eigenfrequencies of the linearized system versus rotor speed σ for n ˜ = 2. Blue: α → 0, Green: α = 0.005, Red: α = 0.01 from a tautochronic path. The definition of the path used here is described in Section 1.1. Moreover, the scaled rotor speed ν is close to unity, motivating the definition of a coordinate that represents the small fluctuations in the mean rotor speed, given by ν(z) = 1 + 2 ρ(z) + O ⇒ ν(z)ν (z) = 2 ρ (z) + O 4 4 (3.14) . (3.15) Using the scaling scheme in Equation (3.13) in Equations (3.10 - 3.12) and solving for ρ in the expansion in Equation (3.15), we obtain an expression which can be inserted into the equations that govern the absorber dynamics in Equation (3.18) and the rotor relative deflection dynamics in Equation (3.19) , resulting in a pair of coupled equations for s and w. These are presented below, following a few more assumptions. Note that the absorber is tuned to address the rotor resonance which occurs in the neighborhood of σ = 1/n for an order n torque, since σ = ωΩ (see Table 3.1). To investigate 0 the behavior of the system near this operating point, which is precisely the eigenfrequency 44 veering zone, two detuning parameters δ1,2 are introduced via n ˜ 2 = n2 + δ1 (3.16) 1 σ2 = 2 , n + δ2 (3.17) which completes the scaling assumptions. The first of these corresponds to the detuning of the absorber from the excitation order and the second is associated with variations in the speed near the resonance. These assumptions result in the model to be considered for analysis, which is a set of two weakly coupled, weakly nonlinear, internally resonant, and resonantly driven oscillators of the form ˆ δ1 ) + O s + n2 s = f1 (s, s , wˆ , µ ˆa , φ, 2 ˆ δ2 , z) + O wˆ + n2 wˆ = f2 (s, s , w, ˆ wˆ , σ 2 , Γ, (3.18) 2 (3.19) where f1 = 2s3 φˆ − sδ1 + n2 g(s)wˆ − µ ˆa s − 2n2 swˆ f2 = 1 2 n2 g(s)s − 2δ2 wˆ2 n2 ss − s 2 dg(s) ˆ − Γ sin(nz) ds (3.20) (3.21) These equations are amenable to perturbation techniques, and here we use the method of averaging. The standard polar form representation for oscillations of s and w are employed, as follows: s(z) = a1 (z) cos(nz + β1 (z)) (3.22) s (z) = −a1 (z)n sin(nz + β1 (z)) (3.23) w(z) ˆ = a2 (z) cos(nz + β2 (z)) (3.24) wˆ (z) = −a2 (z)n sin(nz + β2 (z)), (3.25) which implies the standard constraint equations on s (z) and wˆ (z) [27]. These are substituted into Equations (3.18 - 3.21), and the resulting set of equations are then used to solve 45 for a1,2 and β1,2 . Those results are slowly varying in time with small oscillating terms, are they are averaged over one period of z, that is, over 2π/n, resulting in equations that govern the slowly varying amplitudes and phases. These are given by a ¯1 = − a ¯1 β¯1 = − a ¯2 = − a ¯2 β¯2 = − 2π n n 0 2π n n 0 2π n n 0 2π n n 0 f1 sin(nz + β1 (z))dz (3.26) f1 cos(nz + β1 (z))dz (3.27) f2 sin(nz + β2 (z))dz (3.28) f2 cos(nz + β2 (z))dz. (3.29) Here the integrals are complicated and are not evaluated in closed form. The equilibrium solutions, given by the zeros of these functions, and their stability, are calculated by numerical evaluation of these integrals. This allows one to efficiently carry out parameter studies that include stability information. 3.4 Case Studies The averaged equations of motion of the system, which govern the slowly varying amplitudes and phases, are numerically solved for equilibria in the vicinity of the eigenvalue veering, that is, near resonance. The local stability of these equilibrium points is determined by the eigenvalues of the Jacobian of these equations. We consider a small sampling of response curves that demonstrate the method and point out some general trends. Certain system parameters are kept constant throughout the present case studies, namely the excitation order, n = 2, the absorber damping, µ ˆa = 0.05, and the scaling parameter, which is the ratio of absorber inertia to rotor inertia, = 0.1. Assuming a fixed excitation order is realistic since it is set by the configuration of the rotating machine. For example, in a four-stroke internal combustion engine with N cylinders the dominant order of excitation is n = N/2, although higher order harmonics also exist and are generally important for auto- 46 motive crankshaft resonance considerations. Also, these absorbers have much less rotational inertia than the rotor, and since they stay tuned at all rotor speeds, they are more effective when lightly damped, so these assumptions are similarly valid in applications. Without the addition of a CPVA, for an excitation order n the resonance of the rotor system occurs at a dimensionless speed of σ = 1/n, since σ = ωΩ (see Table 3.1). So, for the 0 parameters of interest, the response analysis is first carried out at the dimensionless rotor ˆ and the linear absorber speed σ = 0.5. At this speed, the excitation torque amplitude, Γ, tuning, n ˜ , are varied for paths that are tautochronic, softening, and hardening, in terms of their nonlinear tuning. The response branches for the amplitudes of the CPVA and the rotor deflection are plotted in Figures 3.4b, 3.5b, and 3.6b). The stability of these branches is depicted using color where blue curves are stable and red curves are unstable. For all three nonlinear path options, it can be seen that near n ˜ = 2, the shaft deflection is significantly reduced without a steep increase in the case of detuning. In the cases of softening and hardening paths, it can be seen that the n ˜ value corresponding to minimum shaft deflection shifts as the torque amplitude in increased, as expected. Moreover, in all three cases, unstable zones are observed to appear and widen as the torque amplitude is increased. At first, it seems logical to off tune the CPVA at the expense of reduced suppression, so that the unstable zones are not within the operation range. This is the most common tuning approach when the CPVAs are used to reduce rotor rigid body torsional vibration. However, one also needs to inspect the behavior at neighboring rotation speeds before selecting the tuning. The next two sets of parameter studies are at slightly slower and slightly higher dimensionless rotor speeds. The excitation torque amplitude steps used above are repeated in this example. The response amplitudes obtained in these cases are plotted in Figures 3.4 - 3.6. The regions of increased CPVA response are observed at tunings slightly higher than n ˜=2 for σ = 0.48 and at slightly lower than n ˜ = 2 for σ = 0.52, when compared to the results from the horizontal eigenvalue branches obtained from the linear eigenvalue analysis. The asymptotes of these branches separate from each other as the inertia ratio, α, increases. 47 Thus, they approach n ˜ = n as α → 0. The shaft deflection amplitude caused by these peaks, on the other hand, is much lower when compared to the reduced shaft resonance response seen in Figures 3.4b, 3.5b, and 3.6b. 48 0.25 0.2 0.2 0.2 0.15 0.15 0.15 a1 0.25 a1 a1 0.25 0.1 0.1 0.1 0.05 0.05 0.05 0 1.5 2 n ˜ 0 1.5 2.5 0.3 2 n ˜ 0 1.5 2.5 2 2 n ˜ 2.5 2 n ˜ 2.5 0.25 0.2 1.5 0.2 a2 a2 a2 0.15 1 0.1 0.1 0.5 0 1.5 2 n ˜ (a) σ = 0.48 2.5 0.05 0 1.5 2 n ˜ (b) σ = 0.50 2.5 0 1.5 (c) σ = 0.52 Figure 3.4: Tautochronic path steady state CPVA (a1 ) and rotor deflection (a2 ) steady state amplitude values. Blue: stable, ˆ = 0.01 to 0.85. Red: unstable, Red-dashed: CPVA cusp amplitude. n = 2, = 0.1, φ = 0, Γ 49 0.25 0.2 0.2 0.2 0.15 0.15 0.15 a1 0.25 a1 a1 0.25 0.1 0.1 0.1 0.05 0.05 0.05 0 1.5 2 n ˜ 0 1.5 2.5 0.3 2 n ˜ 0 1.5 2.5 2 2 n ˜ 2.5 2 n ˜ 2.5 0.25 0.2 1.5 0.2 a2 a2 a2 0.15 1 0.1 0.1 0.5 0 1.5 2 n ˜ (a) σ = 0.48 2.5 0.05 0 1.5 2 n ˜ (b) σ = 0.50 2.5 0 1.5 (c) σ = 0.52 Figure 3.5: Softening path steady state CPVA (a1 ) and rotor deflection (a2 ) steady state amplitude values. Blue: stable, ˆ = 0.01 to 0.85. Red: unstable, Red-dashed: CPVA cusp amplitude. n = 2, = 0.1, φ = −5, Γ 50 0.25 0.2 0.2 0.2 0.15 0.15 0.15 a1 0.25 a1 a1 0.25 0.1 0.1 0.1 0.05 0.05 0.05 0 1.5 2 n ˜ 0 1.5 2.5 0.3 2 n ˜ 0 1.5 2.5 2 2 n ˜ 2.5 2 n ˜ 2.5 0.25 0.2 1.5 0.2 a2 a2 a2 0.15 1 0.1 0.1 0.5 0 1.5 2 n ˜ (a) σ = 0.48 2.5 0.05 0 1.5 2 n ˜ (b) σ = 0.50 2.5 0 1.5 (c) σ = 0.52 Figure 3.6: Hardening path steady state CPVA (a1 ) and rotor deflection (a2 ) steady state amplitude values. Blue: stable, ˆ = 0.01 to 0.85. Red: unstable, Red-dashed: CPVA cusp amplitude. n = 2, = 0.1, φ = 5 , Γ 51 3.5 Conclusion In this chapter, the nonlinear dynamic response of a flexible shaft equipped with a centrifugal pendulum vibration absorber is investigated. Previous work [36] demonstrated the linear eigenvalue structure of the system and showed the possibility of eliminating shaft resonance using a properly tuned CPVA. This requires a special tuning that is based on the eigenvalue veering behavior encountered as the rotation speed is varied. In practice, CPVAs move through large amplitudes and behave nonlinearly, and these effects must be taken into account to implement this tuning practice, and this is the purpose of the analysis presented in this chapter. An idealized model that captures the essential dynamics is considered. The equations of motion of this model are put into a dimensionless form and certain physical quantities are scaled according to expected orders of magnitudes for these quantities. As a result, a set of equations amenable to the method of averaging is obtained. Steady state amplitude and phase conditions are determined from the averaged equations by numerical integration. A set of parameter studies is carried out in the neighborhood of the eigenvalue veering zone, that is, near resonance, which is where the system robustness to nonlinear effects must be considered. It is shown that the response amplitudes have the desired characteristics, but that several instability zones were detected that cannot be predicted by a linear analysis. In previous studies that considered a rigid rotor model, the allowable maximum torque level was dictated by either these instability zones or by kinematic limits of the path. Generally, in those cases the absorber was off-tuned at the expense of reduced performance in order to increase the excitation limit. In contrast, for the present model, which considers absorbers tuned to address a torsional resonance, the upper limit for the level of off-tuning selected is set by the veering branches of the eigenvalues, which, if encountered by a engine order line, would result in a resonance. Increasing the effective inertia of the CPVA can shift this upper limit such that resonance can be avoided, or the system can be designed so 52 that the maximum allowed excitation torque is determined not by the absorber amplitude limits, but by the stability boundaries of the desired response. In addition, one might avoid these instabilities if some absorber damping can be introduced into the system. All these solutions, as usual, will reduce the effectiveness of the absorber. However, depending on the nature of the instability, i.e., whether the bifurcation is sub- or super-critical, one may be able to operate into the unstable zones, at least for short periods of time. These issues, and their relevance to practical designs, are left for future studies. 53 CHAPTER 4 INCREASING THE EFFECTIVE INERTIA OF BIFILAR SUSPENSION CPVAS BY EMPLOYING THEIR ROTATIONAL INERTIA THROUGH KINEMATIC ALTERATIONS 4.1 Introduction (a) (b) (c) Figure 4.1: Illustration of different types of CPVA suspension designs. (a) Ring/Roller Type, (b) Compound Pendulum Type, (c) Bifilar Suspension Type Implementing a centrifugal pendulum vibration absorber requires a suspension mechanism that prescribes the relative motion between the absorber and the host rotor. Many designs have been developed to suspend the pendulum absorber mass to the rotor of interest, including simple pendulums, ring/roller absorbers, and most commonly a bifilar (two point) suspension, which employs rollers moving along pairs of identical, but inverted, tracks cut from the pendulum and the rotor. These suspensions have been known for decades and are described in [10], some of which are illustrated in Figure 4.1. The bifilar mechanism uses two sets of cutouts on the absorber and on the rotor to guide the motion of the absorber along its prescribed path using rollers. This mechanism allows one to use relatively large absorber masses within the confined spaces of engines, and permits the absorber to follow noncircular paths, as required for proper tuning of the nonlinear dynamics of the absorber at large amplitudes [14]. However, the standard bifilar suspension prevents relative rotation 54 between the absorber and the rotor, and thus its rotational inertia about its center of mass (COM) acts as a dead weight that is fixed on the rotor, while its mass acts like a point mass moving along the COM. Figure 4.2: Conceptual drawing of the system. Here the black outlined boxes represent the motion of the absorber for a standard bifilar mechanism while the red dashed boxes represent the positions of the rocking absorber. S is the arclength coordinate that indicates the position of the COM of the absorber, and ϕ(S) is the rocking angle of the absorber as a function of S. In this study, we investigate a new bifilar suspension mechanism design that allows one to prescribe both translation and rotation of the absorber relative to the rotor. The concept of this rocking absorber behavior is depicted in Figure 4.2. The modification included here imposes a certain rotational motion to the absorber by creating the cutouts accordingly. This allows the rotational inertia of the absorber to join the effort in counteracting the excitation, eventually getting more correction for the same amount of mass when the rest of the conditions are the same. As the rate of rotation imposed in this method is controlled by the cutouts, it becomes another design parameter, similar to the linear and nonlinear tuning path parameters, and the absorber mass and placement. It is well known that the effective inertia of a CPVA, that is, the mass times the square of the distance of the absorber COM to the rotation center for a point mass case, has a major role on its performance. The behavior of the CPVA can simply be imagined as a synchronous energy storage device that works out of phase with the forcing acting on the 55 rotor. As the periodic forcing accelerates the rotor, the CPVA absorbs the energy partially, thus effectively counteracting the forcing. As the forcing tries to decelerate the rotor, the stored energy is fed back to rotor. This out of phase behavior relies on the fact that the order of periodic forcing matches the tuning of the absorber. In turn, in any absorber tuning problem, available absorber inertia plays an important role; the more inertia the absorber has, the more energy it can absorb and return as it moves. Generally, absorbers are detuned from the excitation order as a trade-off in performance, so that they can handle a large range of excitation amplitudes, and so that they are robust to imperfections [7, 16]. The selection of this detuning is highly correlated to how much absorber inertia is available [7]. In most applications, the available space is a key limiting factor on how much absorber inertia can be included in a given design. Another consideration is the amount of mass and rotational inertia added to the rotor by the absorbers, which affects system responsiveness, weight, and cost. Employing the rotational inertia of the absorber in order to increase the effective inertia, without adding any extra mass, is therefore highly desirable in practical applications. This design was independently conceived by the author, although similar technologies had very recently been proposed in the literature [37], and had been numerically analyzed in [38]. In fact, such designs are a natural blend of absorbers that roll and those that translate, for example, as shown in Figure 4.1. In the present work we provide the first systematic, quantitative analysis of these systems based on their kinematics, and provide a detailed analysis of their nonlinear response, which is required to fully understand and appreciate the benefits of these designs. In the following sections, we elaborate on the construction of the cutouts that guide the bifilar suspension rollers for these so-called “rocking” CPVAs. Through these cutouts it becomes possible to independently prescribe the COM translation path and the constrained rocking motion of the absorber mass. We next discuss the tuning of these absorbers and present the dynamical analysis framework for this designs. Finally, we show analysis results 56 for a number of sample system configurations that reveal the effects of rocking motion on the performance of CPVAs, and compare them with their purely translational counterparts. 4.2 Cutout Kinematics As stated previously, the typical bifilar suspension mechanism employs two rollers between the absorber and the rotor. Even tough these rollers are free to roll along the cutout surfaces, in operation they always stay in contact with the absorber and the rotor and no extra degrees of freedom are created, provided that the tangents of the contact surfaces are always perpendicular to the normal force due to the centrifugal field and that this condition ensures that the rollers roll without slip relative to their contact surfaces. The cutout formulation for a given COM path is obtained in [14] for a standard bifilar suspension. In that work, the path followed by the center of a roller is taken to be one half of the path followed by the absorber COM. The cutouts are then obtained by offsetting the roller path by its radius. The principle of this approach can be illustrated for a flat object moving parallel to another flat surface with a roller between them. If the motion is rolling without slip, the translation of the roller center is one half that of the translation of the flat object, as shown in Figure 4.3 . The reason this principle can directly be carried over to the standard bifilar mechanism is that the only relative motion between the absorber and the rotor is planar translation. ∆x ∆x/2 Figure 4.3: Illustration of moving a mass over a flat surface with a roller. 57 However, the situation is different for the proposed version of the bifilar mechanism. Contrary to the translation-only design, it is desired to impose a rotation of the absorber about its COM, resulting in a motion in which the direction and amplitude of the velocity vectors on the absorber depend on the local coordinates of the points on the absorber, as shown in Figure 4.4 . θ S Figure 4.4: Vectors indicating that each point on the absorber will have different velocities as a result of the combined translation of the COM and rotation about the COM. This implies that, as the motion of the absorber progresses along its COM path, the roller contact point will have a specific direction to move that depends on the history of the rotation moving away from the vertex position. For a given COM path and rotation rate, the cutouts are obtained by a recursive algorithm, with which for every step, the local coordinates of the roller contact point are calculated and used as the seed for the direction and magnitude of motion of that point. The roller center is then moved in that direction by one half of the progression of the absorber at that point. The corresponding cutout points are also marked at a distance equal to the roller center in a direction perpendicular to the direction of motion. In the present work, we consider the case where the absorber rotation is a cubic polynomial in the absorber translation distance S, although this relationship can be easily generalized. 58 In Figure 4.5 we show examples of the cutouts obtained using this approach for the same COM path and the same vertex roller positions, for four values of the rotation rate. In the top-left case, the cutouts are calculated for a non-rocking bifilar mechanism, and in the bottom-right case the rocking rate is chosen so that the rotation at the cusp point is 21◦ . We now turn to an analysis of the dynamics of a rotor-CPVA system with these kinematics. -0.11 -0.11 -0.12 -0.12 y -0.1 y -0.1 -0.13 -0.13 -0.14 -0.14 -0.15 -0.04 -0.02 0 x 0.02 -0.15 -0.04 0.04 -0.11 -0.11 -0.12 -0.12 0 x 0.02 0.04 -0.02 0 x 0.02 0.04 y -0.1 y -0.1 -0.02 -0.13 -0.13 -0.14 -0.14 -0.15 -0.04 -0.02 0 x 0.02 0.04 -0.15 -0.04 Figure 4.5: Samples of cutouts for rocking bifilar absorbers for the same COM path with different rocking rates. Here, the red and green dashed lines represent cutouts on the rotor while red and green solid lines represent the absorber cutouts. The solid outlined circles represent the rollers shown at their center position, and the blue dashed line is the COM path. The rocking rates, moving from top left to bottom right, are set such that the rocking angle at the cusp becomes 0◦ , 7◦ , 14◦ , 21◦ . 59 4.3 Dynamical Model In this section we describe results from analysis of the system response using both the method of averaging and a generalization of the extended harmonic balance method described in Chapter (2). Perturbation methods are the preferred approach in the study of dynamics of pendulum absorbers, since they produce analytical results that provide explicit parameter dependence as well as physical insight. As long as the small parameter selection is conducted in a way that the nonlinearities in the equations of motion are represented accurately as small perturbations, these techniques allow for an accurate assessment of the system dynamics. On the other hand, the HBM method requires tedious variable recasting during the initial setup of the analysis, and is purely numerical in nature, but it is more generally applicable to a wider range of systems. However, once the system of equations are reshaped into the required form, it can be used to obtain the dynamic response as accurately as numerical integration methods, while being able to do so as efficiently as perturbation methods. The first step to either of these solution methods, as well as to the crucial absorber tuning analysis, is the derivation of the system equations of motion. The system investigated here is modeled using Lagrange’s equations. The Lagrangian can be expressed accurately by the total kinetic energy of the rotor-CPVA system, since the rotation of the system dominates the effects of gravity. The Lagrangian of the system reads as follows: 1 1 L = J θ˙2 + mj 2 2 N 2 Rj2 (Sj )θ˙2 + S˙ j2 + 2Gj (Sj )θ˙S˙ j + rg,j θ˙ + ϕj (Sj )S˙ j 2 (4.1) j=1 Here, θ is the angular position of the rotor and Sj is the (arc length) displacement of the j th absorber COM relative to the vertex on its path. Parameters mj and rg,j are the mass and radius of gyration of absorber j, respectively. The prescribed rocking angle of the absorber is given by ϕj (Sj ), which can be quite general. In this study we take ϕj to be a cubic polynomial in Sj with only odd powers, in order to obtain a symmetric rocking motion 60 about the vertex, as follows, 1 ϕj (Sj ) = B0,j Sj + B1,j Sj3 . 3 (4.2) Terms related to functions Rj (Sj ) and Gj (Sj ) are related to the COM path, as explained in the Introduction chapter. Note that the kinetic energy of the system is modified by the rocking absorber motion, with the addition of only a few terms that appear at the end of the energy expression. The Lagrangian equations of motion of this system, with the applied torque and damping forces added, read as follows,   N J + mj 2 Rj2 (Sj ) + rg,j 2 B0,j + B1,j Sj2 mj Gj (Sj ) + rg,j S¨j j=1 j=1 dRj2 (Sj ) N mj + N  θ¨ + dSj j=1 dGj (Sj ) 2 2 S S˙ 2 θ˙S˙ j + S˙ j + 2B1,j rg,j j j dSj N = T0 + T (t) − c0 θ˙ + ca,j Gj (Sj )S˙ j (4.3) j=1 2 mj Gj (Sj ) + rg,j B0,j + B1,j Sj2 1 − mj 2 2 θ¨ + mj 1 + rg,j B0,j + B1,j Sj2 2 S¨j dRj2 (Sj ) 2 2 ˙2 2 B θ˙ + 4rg,j 1,j B0,j + B1,j Sj Sj dSj = −ca,j S˙ j (j = 1 to N ) (4.4) where T (θ) is the applied torque acting on the rotor, and ca,j represents the equivalent viscous damping coefficient for absorber j as it moves along its path, T0 is the constant torque applied to the rotor in order to maintain a mean operating rotational speed Ω, by counteracting viscous rotor damping that is modeled by coefficient c0 . The forcing term T (θ) in Equation (4.3) has the form of engine order excitation, so it is most naturally expressed as a function of rotor position, and it is also assumed to be dominated by a single harmonic, so that we take T (t) = T sin(nθ) 61 (4.5) where n is the engine excitation order. 4.3.1 Non-dimensionalization Table 4.1: Definition of dimensionless variables. Variables αj Definitions mj Rj2 (0)/J Variables βj Definitions rg,j/R (0) j b0,j B0,j Rj (0) b1,j B1,j Rj3 (0) µa,j ca,j/m Ω j sj Sj/R (0) j gj (sj ) Gj (Sj )/R (0) j rj2 (sj ) Rj2 (Sj )/R2 (0) µ0 c0/JΩ Γ0 T0/JΩ2 Γ T /JΩ2 j The equations of motion are non-dimensionalized using the definitions listed in Table (4.1). As none of the terms explicitly depend on time in the equations, it is convenient to replace time as the independent variable with the rotor angular position θ. Also, it is convenient to introduce the following definition for the dimensionless mean speed of the rotor: ν(θ) = θ˙ . Ω (4.6) These operations and definitions, along with use of the chain rule, yield the following dimensionless set of equations that govern the behavior of the rotor and the absorbers: 62   N αj βj2 + rj2 (sj ) + gj (sj ) + βj2 b0,j + b1,j s2j 1 + sj  νν j=1   N αj gj (sj ) + βj2 b0,j + b1,j s2j +  ν 2s j j=1  N + αj j=1 drj2 (sj ) dsj  dgj (sj ) + 2b1,j βj2 sj sj  ν 2 sj dsj   + N = Γ0 + Γ(θ) − µ0 ν +  αj µa,j gj (sj )sj  ν (4.7) j=1 gj (sj ) + sj + βj2 b0,j + b1,j s2j + 1 + b0,j + b1,j s2j 2 1 + b0,j + b1,j s2j sj ν νsj 2 1 drj (sj ) − ν + 2βj2 b1,j b0,j + b1,j s2j νsj sj = −µa,j sj 2 dsj (j = 1 to N ) (4.8) d() where () = dθ . 4.3.2 Absorber Tuning In the case of a classical, non-rocking bifilar suspension implementation, which causes the absorbers to act as point masses, the absorber tuning, both linear and nonlinear, is controlled solely by the parameters that dictate the path followed by the absorber COM. That is, the tuning is determined by considering constant rotor speed, and the absorber inertial properties vanish from the absorber equation of motion in this case. r2 (s) = 1 + φ0 s2 + φ1 s4 , 63 (4.9) where φ0 = 1 − ρˆ−1 0 φ1 = − (4.10) −1 + λ2 + ρˆ0 12ˆ ρ30 with ρˆ0 := ρ0/R(0). (4.11) As shown in the path definition, Equation (4.9), the coefficient of the quadratic term, φ0 , in r2 (s) depends solely on the radius of curvature, ρ0 , of the path at the vertex. The quartic term, φ1 , in r2 (s), on the other hand, depends on both ρ0 and the nonlinear path parameter λ. Since this higher order term in the absorber path vanishes in the linearized dynamics, only φ0 controls the small-amplitude (linear) tuning. So, one determines ρ0 to obtain a desired linear tuning and then adjusts λ to set the amplitude dependent (nonlinear) tuning, including the special case of tautochronic (amplitude independent) tuning. To see a more detailed definition on paths, see Section 1.1. In the case of a rocking CPVA, the tuning depends not only on the path parameters, but also on the rocking parameters and the radius of gyration of the CPVA mass. Here we derive the expressions that yield the amplitude dependent tuning of these rocking CPVAs using a model consisting of a rotor and a single rocking absorber. First, we write the Hamiltonian of a single undamped absorber at constant rotor speed, as HΩ = S˙ ∂LΩ − LΩ ∂ S˙ (4.12) where LΩ = L 1 1 = JΩ2 + m 2 2 ˙ θ=Ω R2 (S)Ω2 + S˙ 2 + 2G(S)ΩS˙ + rg2 Ω + ϕ (S)S˙ 2 . (4.13) Applying the dimensionless parameter substitutions listed in Table (4.1), we obtain 2 1 HΩ = αJ −Ω2 β 2 + r2 (s) + 1 + β 2 b0 + b1 s2 s˙ 2 . 2 (4.14) Note that generally the Hamiltonian is expressed in terms of the generalized coordinates and momenta, but here it is convenient to keep s˙ instead of replacing it with the associated momentum variable. 64 This Hamiltonian formulation approach allows one to extract the relation between the instantaneous absorber velocity and the position of the CPVA under the assumption of constant rotor speed [39]. This can then be integrated to determine the absorber period as a function of amplitude. Note that, this single degree of freedom system is conservative. Hence, H˙ Ω = 0 or HΩ = const., which implies HΩ = HΩ = const. (4.15) s=A,s=0 ˙ Therefore, ds =Ω dt =Ω r2 (s) − r2 (A) 1 + β 2 b0 + b1 s 2 (4.16) 2 φ0 (s2 − A2 ) + φ1 (s4 − A4 ) (4.17) 2 1 + β 2 b0 + b1 s 2 From this, we can find the period of oscillation as a function of the amplitude of oscillation, A, by calculating the following integral,   2   A 2 2 1 + β b 0 + b1 s 4 T (A) = ds Ω 0  φ0 (s2 − A2 ) + φ1 (s4 − A4 )  (4.18) Numerically evaluating this integral using a trapezoid scheme presents a difficulty, as the integrand has a pole at the upper integration limit, s = A. Moreover, expanding the integrand in a series summation also fails to accurately approximate the integrand near the integration limit. To overcome these integration problems, first we expand only the numerator of the integrand, which is always positive, in a series summation in s up to order M , as 1 + β2 b 0 + b1 2 s2 M ak s k . = (4.19) k=0 Then, the integral takes the following form, which can be rewritten as combinations of complete elliptic integrals of the first and second kind [40]. Then, the integral can be evaluated robustly, that is, without difficulties at the upper limit. 4 A T (A) = Ω 0 M k k=0 ak s φ0 (s2 − A2 ) + φ1 (s4 65 − A4 ) ds. (4.20) Note that the result is linear in the rotor mean speed, Ω, for any rocking function, indicating that order tuning is valid for a rocking CPVA. As per its definition, the absorber order tuning is obtained by taking the reciprocal of the absorber period and dividing it by Ω, as n ˜ (A) = 2π . ΩT (A) (4.21) It has been found that the integral in Eqn. (4.20) yields quite different expressions in cases with nonlinear parameter φ1 = 0, which yields a tautochronic tuning for a non-rocking CPVA, when compared to those with φ1 = 0. Even though both cases result in closed form expressions, the latter involves complete elliptic integrals of the first and second kinds. Tuning of the absorbers in both cases are given below for M = 4, which is sufficiently accurate for most applications, 1 + β2 b0 + b1 2 s2 M ak s k = (4.22) k=0 ≈ a0 + a2 s2 + a4 s4 (4.23) where a0 = 1 + β 2 b20 1/2 a2 = β 2 b0 b1 1 + β 2 b20 (4.24) −1/2 −3/2 1 a4 = β 2 b21 1 + β 2 b20 . 2 (4.25) (4.26) Substitution of the expansion Equation (4.22) into Equation (4.21) yields the following amplitude dependent tuning expression, √ 3πφ21 −F3 n ˜ (A) = 2 F1 −3a0 φ21 − F3 2a4 φ0 − 3a2 φ1 + A2 a4 φ1 + F2 F3 (2a4 φ0 − 3a2 φ1 ) (4.27) where F1 = K F2 = E φ0 −1 A2 φ1 + φ0 φ0 −1 2 A φ1 + φ0 F3 = φ0 + A2 φ1 66 (4.28) (4.29) (4.30) with K(x) and E(x) being the complete elliptic integrals of the first and second kinds, respectively, with the elliptic modulus of x. For the case φ1 = 0, taking the COM path to be tautochronic in the non-rocking case, this result simplifies significantly to √ 8 −φ0 n ˜ (A) = 8a0 + 4A2 a2 + 3A4 a4 4.3.2.1 (4.31) Linear Tuning Small amplitude tuning is a primary metric in the design of CPVAs. Here, using the expressions obtained by the integration of the absorber period, we can extract the small amplitude tuning by just evaluating these at A → 0, yielding n ˜ (0) = 2π ΩT (0) (4.32) = −φ0 1 + β 2 b20 (4.33) = 1 − ρˆ0 . ρˆ0 1 + β 2 b20 (4.34) It turns out both Equation (4.27) or Equation (4.31) yield the same result for the linear tuning, which is expected since φ1 is a nonlinear path parameter that does not affect the small amplitude behavior of the absorber. By rearranging Equation (4.32), one can extract the vertex radius of curvature needed for a given linear tuning order n ˜ (0) as ρ0 = R(0) 1 + n ˜ (0)2 1 + β 2 b20 −1 , (4.35) which is seen to depend on the linear rocking rate, b0 , and the absorber radius of gyration, β. However, the cubic rocking rate does not have any effect on the linear tuning, as expected. This gives the designer the liberty of tailoring the rocking motion such that a given amount of translation and rocking space for an absorber can be utilized as needed by the application at hand. First, Equation (4.32) suggests that, for a fixed translation path, an increase of the linear rocking rate causes the linear tuning order to reduce. Thus, one needs 67 to adjust the translation path accordingly. This, however, would also cause the cusp point of the path to change. On the other hand, the designer might choose to utilize the entire rocking space available with a purely cubic rocking rate, which implies that the performance gain from rocking is induced only at higher absorber amplitudes. In essence, having these individual parameters at the disposal of the designer widens the design space for a given set of constraints. In order to reveal the trade-offs of these parameters, a detailed dynamical response analysis needs to be conducted in conjunction with the tuning studies. 4.3.2.2 Tautochronic Tuning Tautochronic (amplitude independent) tuning of CPVA systems is a preferred starting point in their design, since it avoids nonlinear effects to leading order. In addition, choosing a near or exact tautochronic path simplifies the perturbation analysis for dynamical response characteristics. With the presence of rocking behavior, one can achieve exact tautochronic tuning only if the rocking rate of the absorber is linearly dependent on the COM path coordinate s. This can be seen in Equation (4.31), where the only way to eliminate the absorber amplitude, A, from the tuning is to set the cubic rocking rate, b1 , to zero. That is, for β = 0, (˜ n(A) = n ˜ (0)) ⇔ (a2 = 0 ∧ a4 = 0) ⇔ (b1 = 0) (4.36) where ∧ implies that both conditions must hold. In the more general case with φ1 = 0, the amplitude dependent tuning expression is given by Equation (4.27). As the term F3 in the numerator of Equation (4.27) is a function of A unless φ1 = 0, one can not achieve amplitude independent tuning with rocking absorbers under the condition of φ1 = 0. 68 4.3.3 Scaling and Averaging A scaling scheme is applied to the equations of motion so that the contribution of terms that are small in practice fall at the correct order for the purpose of applying a perturbation analysis. To this end, the following scaling definitions are used in Equations (4.7 and 4.8), where , the ratio of absorber inertia to rotor inertia, is the central small parameter, αj = α ˆ j , βj = √ ˆ µa,j = µ βˆj , Γ = Γ, ˆa,j (4.37) which captures the small nature of the dimensionless parameters associated with the absorber inertia (both mass and rotational inertia) relative to the rotor inertia, the absorber damping, and the applied torque. Moreover, the scaled rotor speed ν is close to unity, motivating the definition of a coordinate that represents the small fluctuations in the mean rotor speed, given by ν(θ) = 1 + ρ(θ) + O ⇒ ν(z)ν (θ) = ρ (θ) + O 2 2 (4.38) . (4.39) Using the scaling scheme in Equation (4.37) in Equations (4.7 and 4.8), and solving for ρ in the expansion in Equation (4.39), we obtain an expression which can be inserted into the equations that govern the absorber dynamics. This is the uncoupling of the rotor dynamics from the absorber dynamics, as done in many previous studies of these systems. In addition, two detuning parameters, δj,1 and δj,2 are introduced, as follows, n ˜ 2j = n2 + δj,1 φ1,j = δj,2 , (4.40) (4.41) which account for the facts that the tuning of the absorber is close to the excitation order n, and that the nonlinear absorber tuning is close to tautochronic. Using these definitions for the tuning yields the following path, as described in the linear tuning discussion above, r12 (s1 ) = 1 − 1 + b20 βˆj2 n2 + δj,1 s2 + δj,2 s4 . 69 (4.42) This completes the scaling assumptions. These assumptions and expansions in result in the model to be considered for analysis in the single absorber case, i.e. j = 1, which is a forced oscillator with weak nonlinearity, weak near-resonant excitation, and weak damping, of the form s1 + n2 s1 = f1 + O 2 (4.43) where ˆ sin nθ −g1 (s1 ) − s − µ ˆa,1 s1 − δ1,1 s1 + 2δ1,2 s31 f1 = Γ 1 − 2βˆ12 b21,1 s31 s12 − 2βˆ12 b0,1 b1,1 s1 s12 + α ˆ 1 s13 g1 (s1 ) + α ˆ 1 g1 (s1 )g1 (s1 )s12 + n2 s1 βˆ12 b1,1 s1 2 b1,1 s1 2 + 2b0,1 − α ˆ 1 3g1 (s1 ) s1 + g1 (s1 ) 2 + 2s1 2 . (4.44) We can drop the absorber index for clarity, as the remainder of the analysis is conducted for system with a single absorber. This equation is amenable to perturbation techniques, and here we use the method of averaging. The standard polar form representation for the driven response of s is employed, as follows: s(θ) = a(θ) cos(nθ + ψ(θ)) (4.45) s (θ) = −a(θ)n sin(nθ + ψ(θ)), (4.46) which imposes the standard constraint equation on s (θ) [27]. These expressions are substituted into Equation (4.43) , and the resulting set of equations are then used to solve for (a, ψ). Those results are then averaged over one period of θ, that is, over 2π/n, resulting in equations that govern the slowly varying amplitudes and phases, which are given by a ¯ =− a ¯ψ¯ = − 2π n n 0 2π n n 0 f sin(nθ + ψ(θ))dθ (4.47) f cos(nθ + ψ(θ))dθ (4.48) 70 The equilibrium solutions, given by the zeros of these functions, are calculated by numerical evaluation of these integrals. This allows one to efficiently carry out parameter studies, including stability considerations. The averaged equations of motion of the system, which govern the slowly varying amplitudes and phases, are numerically solved for equilibria as a function of system and excitation parameters in order to determine the performance of the system. The steady state solutions, i.e., the zeros of Equations (4.47 and 4.48), for the values of a ¯ and ψ¯ are substituted into the transformation in Equation (4.45) in order to obtain the steady state response of the absorber. This response can then be used in Equation (4.7), along with the expansion in Equation (4.39), to solve for ρ (θ), an approximate representation of the rotor response. This allows one to quantify the reduction in rotor torsional vibration, hence to determine the performance of the absorber. Since the absorber is designed to operate on order n harmonics, we pass ρ (θ) through a Fourier integral to obtain the harmonic n magnitude of the rotor response by ρn = n π 2π/n e−inθ ρ (θ)dθ . (4.49) 0 We now turn to some case studies that demonstrate the analysis tools and provide results about the benefits of rocking absorbers. 4.4 Case Studies The interconnected nature of the parameters that control the behavior of the rocking absorbers expands the size of the parameter space considerably. Here, in order to obtain sample results, we select certain parameters to be fixed and investigate the effects of the remaining parameters. We selected fixed values for the dimensionless absorber inertia ratio (α = 0.1) and for the dimensionless radius of gyration (β = 0.5), and for the desired linear tuning order (˜ n(0) = 2.05). We also fixed the value of the dimensionless quartic coefficient for the absorber 71 path (φ1 ). This parameter is generally referred as the nonlinear tuning parameter in CPVA analysis. However, in this setting, the rocking parameters also affect the nonlinear tuning characteristics, and so we avoid referring it as the sole nonlinear tuning parameter. Finally, the excitation order n is also chosen to be fixed as 2.0. Designing for a desired dynamical behavior of the rocking absorber is an iterative procedure, as it is in the classical non-rocking absorber case. However, in the classical case, the number of parameters at the disposal of the designer are small compared to the present system. The basic rules of thumb for tuning non-rocking absorbers specify that the absorbers should be slightly over-tuned in linear range with respect to the excitation order, n ˜ (0) > n, and the nonlinear tuning of the path should be then adjusted to obtain a tautochronic, or a slightly hardening, path, such that potential undesirable nonlinear affects are avoided, and the operating torque range is extended. However, increasing the level of either the linear or nonlinear over-tuning decreases the torsional vibration reduction performance of the absorber [7, 38]. Introduction of the rocking motion to the design of the absorber affects the above mentioned dynamical characteristics for a given choice of path. As discussed in Section 4.3.2, the linear rocking rate, b0 , directly alters the resulting linear tuning. In other words, for a given linear tuning, the required vertex radius of curvature will be different than that for the corresponding non-rocking case, as shown in Equation (4.32). Thus, the procedure of determining the COM path should include consideration of the rocking behavior, since a given path designed for a non-rocking absorber may result in a relatively conservative or aggressive absorber behavior when accompanied with rocking. Here we present the dynamical response characteristics of a small sample of many possible tuning configurations that involve rocking. Specifically, in these cases we varied the linear and cubic rocking rates, b0 and b1 , respectively. In Table 4.2, the cases analyzed are listed. As mentioned, n ˜ (0) and φ1 are held constant. Consequently, the vertex radius of curvature, the parameter λ, and the cusp of the path are determined by these fixed parameters for 72 Table 4.2: Parameter sets for the rocking rate study. Case A B C D E b0 b1 0.00 0.30 0.60 0.30 0.30 0.00 0.00 0.00 -5.00 5.00 λ 0.918 0.919 0.924 0.919 0.919 φ1 scusp -0.40 -0.40 -0.40 -0.40 -0.40 2.12e-01 2.04e-01 1.84e-01 2.04e-01 2.04e-01 each case. As expected, only a change in the linear rocking rate causes any change in these parameters. The resulting rocking angles and amplitude dependent (nonlinear) tuning of the absorbers considered in cases A through E are plotted in Figure 4.6. Note that the angle of rocking is already a dimensionless parameter. Therefore, irrespective of the dimensionality of the parameters, the resulting rocking angle is absolute (see Table 4.1 for the relations between the dimensional and non-dimensional parameters). Also note that, in Figure 4.6 the absorber displacement is shown as normalized with respect to its cusp amplitude. The corresponding cusp values are provided in Table 4.2. The response characteristics are obtained by varying the dimensionless oscillating torque amplitude, Γ, over a range where the corresponding absorber amplitude is below the kinematic singularity of the path, scusp . Here we investigate the effects of the rocking rates, b0 and b1 , on the system dynamics. Order n harmonic amplitudes of the normalized absorber amplitude, s/scusp , and the corresponding rotor acceleration, νν , are plotted versus Γ in Figure 4.7 for the cases displayed in Table 4.2. As can be seen in Figure 4.7, increasing the linear rocking rate, b0 , provides more correction on the rotor torsional vibration. However, this improvement comes at the price of a reduced range of operating torque amplitude. That is, with more rocking rate, the absorber works more aggressively and therefore has a sharper rate of change in amplitude with respect to the excitation torque, and it reaches its maximum sooner. Moreover, as discussed above, an increased linear rocking rate also reduces the path cusp amplitude, and thus the absorber 73 also has less space to move. Moreover, in cases where the cubic rocking rate was varied, it showed that the inclusion of cubic rocking can partially mitigate the shortcomings of a linear-only rocking motion. We used the value of the linear rocking rate from case B and populated cases D and E with added cubic rocking. It is observed that the increase in percent correction from A to B came with a reduction in the maximum allowable Γ. However, in case D, where the rocking angle has an additional negative cubic term, not only is there more correction when compared to the non-rocking case A, but the allowable range of Γ is also increased. This indicates that with the proper adjustment of rocking parameters, the torsional vibration reduction of the absorber can be improved without adding any mass and without sacrificing other performance metrics, such as the maximum allowable excitation amplitude. As a final observation in the study of cubic rocking rates, in this particular setting, we do not see any particular advantage in using a positive cubic rocking rate, such as that considered in case E. 74 8 6 4 ϕ (deg) 2 0 -2 A B C D E -4 -6 -8 -1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 Normalized Absorber Displacement (|s/smax |) Normalized Absorber Displacement (|s/smax |) (a) Rocking angle vs. absorber displacement 1 0.8 0.6 0.4 0.2 0 2.01 2.02 2.03 2.04 2.05 2.06 2.07 2.08 n ˜ (s/smax ) (b) Tuning order vs.absorber amplitude Figure 4.6: Rocking and tuning values of the CPVA for the cases in Table 4.2. 75 Rotor Acceleration Mag. (|ν ν'|) 0.2 Locked Abs. Order n Case A Case B Case C Case D Case E 0.15 0.1 0.05 0 0 0.05 0.1 0.15 0.2 0.25 Excitation Torque (Γ) Absorber Response Mag. (|s/smax|) (a) Rotor acceleration 1 0.8 0.6 0.4 Case A Case B Case C Case D Case E 0.2 0 0 0.05 0.1 0.15 0.2 0.25 Excitation Torque (Γ) (b) Absorber amplitude 90 Case A Case B Case C Case D Case E Percent Correction 80 70 60 50 40 30 0 0.05 0.1 0.15 0.2 0.25 Excitation Torque (Γ) (c) Correction Figure 4.7: Harmonic amplitudes of order n of the response of the systems in the cases from Table 4.2. The black curve in (a) indicates the reference rotor response at order n, when the absorber is locked at its vertex. 76 4.5 Conclusion The present work is based on a model of a rotor equipped with bifilar CPVAs that are kinematically designed to have an imposed rotational motion about their COM that is linked to the translation of the COM. This configuration increases the effective inertia of the absorber, thereby providing better vibration reduction for a given amount of absorber mass. The roller paths and cutouts that will prescribe the desired translational and rotational motions of the absorber are generated numerically, demonstrating how such absorbers can be realized. The governing equations of motion of the corresponding model are put into a dimensionless form and certain physical quantities are scaled according to expected orders of magnitudes for these quantities. As a result, an equation that governs the absorber motion, and is amenable to the method of averaging, is obtained. Steady state amplitude and phase conditions are determined from the averaged equations by numerical methods. A set of parameter studies is carried out to investigate the performance gains achieved by the addition of absorber rotation. The results show that, without adding any extra mass, the reduction in torsional vibrations can be improved by up to 15%. Moreover, it has been found that by using only a linear rocking rate, torsional vibration reduction is improved, but the maximum torsional excitation amplitude at which the absorber can work is reduced, just as in the case of linear over-tuning in the non-rocking case. However, we also found that utilizing a negative cubic rocking rate in conjunction with a positive linear rate, improved the correction performance without sacrificing operating range. This design has the potential for useful practical implementation. This approach includes the linear and cubic rocking rates of the absorber as additional design variables, which allows more flexibility, but more complexity, to the design process. In a given application, parameter studies are needed to find the optimal absorber parameters, which must include both performance and stability considerations. The methods demonstrated in this chapter provide the tools required for such investigations. An extension of 77 the current study would be to consider the stability of the synchronous response of a system of multiple, identical rocking absorbers, similar to that described in [7]. 78 CHAPTER 5 CONCLUSIONS In this dissertation, we contributed to three specific topics in the analysis and design of CPVAs. First, we developed a new framework, based on a harmonic balance method (HBM), for the numerical estimation for the steady state response and stability characteristics of systems equipped with CPVAs. We also posed two new design problems regarding CPVA systems, specifically, their tuning to avoid torsional resonances of the support rotor and increasing absorber inertia by allowing for absorber rotation by so-called “rocking”, and presented detailed modeling and dynamical analysis for these two topics. The nonlinear equations of motion for these systems are generally analyzed either with direct time integration or using perturbation methods, both of which have limitations in terms of speed and/or applicability. The HBM method deveoloped herein is adopted from [21], in which the nonlinearities are put into a quadratic polynomial form. The solution is then expressed as a Fourier series for application of the HBM with a relatively high number of assumed harmonics. Essentially, the method differs from the classical harmonic balance method in that it has a pre-processor for the problem formulation. It operates on the system equations of motion and rewrites them as a larger, yet equivalent, set of equations using auxiliary variables. In this way, the application of the HBM, namely assuming a periodic solution with a given number of harmonics and constructing an algebraic system of equations that govern the amplitudes of the harmonics of the dependent variables, can be automated irrespective of the structure of the original equations of motion. Estimation of the stability of the periodic steady state solutions obtained in this way was also made possible using a related algorithm outlined in [26]. The requirement for this secondary step in the solution process is that each component of the Jacobian of the original equations of motion must be written in terms of its harmonic coefficients at each periodic solution, where 79 stability information is sought. The form of the Jacobian of systems with CPVAs is much more complex than the equations of motion themselves. Therefore, here we evaluated it numerically by reconstructing the periodic solutions of the states of the equations of motion using the coefficients obtained by HBM and extracted the Jacobian harmonic coefficients numerically using a Fast Fourier Transform. This solution framework was tested for two different system configurations that are known to show significant qualitative changes in their response characteristic as the excitation amplitude is varied. The results showed that this approach is efficient and accurate for use in the analysis of CPVA systems. The complexity of the CPVA systems investigated here was kept at a level that makes the demonstration of the method as clear as possible. However, more involved configurations can and have been successfully analyzed with this algorithm. This tool will allow for the rapid study of absorber systems with parameter values that render analytical tools inapplicable and/or direct numerical simulations too slow. Such a tool will be very useful for parameter studies and design purposes. We also considered the nonlinear dynamic response of a flexible shaft equipped with a centrifugal pendulum vibration absorber. Previous work by the author had demonstrated the linear eigenvalue structure of this system and showed the possibility of eliminating shaft resonance using a properly tuned CPVA. This requires a special tuning that is based on the eigenvalue veering behavior encountered as the rotation speed is varied. In practice, CPVAs move through large amplitudes and behave nonlinearly, and these effects must be taken into account to implement this tuning in practice, which motivated the current study. An idealized model that captures the essential dynamics was considered, namely a torsional resonance and a CPVA. The equations of motion of this model are put into a dimensionless form and certain physical quantities are scaled according to expected orders of magnitudes for these quantities. As a result, a set of equations amenable to the method of averaging is obtained. Steady state amplitude and phase conditions were determined from the averaged equations by numerical integration. A set of parameter studies was carried out in the neigh- 80 borhood of the eigenvalue veering zone, that is, near resonance, which is where the system robustness to nonlinear effects must be considered. It is shown that the response amplitudes have the desired characteristics, but that several instability zones were detected that cannot be predicted by a linear analysis. In previous studies that considered a rigid rotor model, the allowable maximum torque level was dictated by either these instability zones or by kinematic limits of the path. Generally, in those cases the absorber was off-tuned at the expense of reduced performance in order to increase the excitation limit. In contrast, for the present model, which considers absorbers tuned to address a torsional resonance, the upper limit for the level of off-tuning selected is set by the veering branches of the eigenvalues, which, if encountered by a engine order line, would result in a resonance. Increasing the effective inertia of the CPVA can shift this upper limit such that resonance can be avoided, or the system can be designed so that the maximum allowed excitation torque is determined not by the absorber amplitude limits, but by the stability boundaries of the desired response. In addition, one might avoid these instabilities if some absorber damping can be introduced into the system. All these solutions, as usual, will reduce the effectiveness of the absorber. However, depending on the nature of the instability, i.e., whether the bifurcation is sub- or super-critical, one may be able to operate into the unstable zones, at least for short periods of time. These issues, and their relevance to practical designs, are left for future studies. We also considered a means of increasing the effective inertia of the absorber mass by allowing it to rotate in a prescribed manner as it moves. To investigate this approach, we developed a model of a rotor equipped with bifilar CPVAs that are kinematically designed to have an imposed rotational motion about their COM that is linked to the translation of the COM. This configuration increases the effective inertia of the absorber, thereby providing better vibration reduction for a given amount of absorber mass. The roller paths and cutouts that will prescribe the desired translational and rotational motions of the absorber are generated numerically. The equations of motion of the model are put into a dimensionless form and certain physical quantities are scaled according to expected orders of magnitudes for 81 these quantities. As a result, an equation that governs the absorber motion, and is amenable to the method of averaging, is obtained. Steady state amplitude and phase conditions are determined from the averaged equations by numerical methods. A set of parameter studies was then carried out to investigate the performance gains achieved by the addition of absorber rotation. The results show that, without adding any extra mass, the reduction in torsional vibrations can be improved by up to 15%, although the torque range over which the system can operate is slightly reduced. In summary, these effects are the result of the fact that the rocking forces the absorbers to work more aggressively. This approach includes the linear and cubic rocking rates of the absorber as additional design variables, which allows more flexibility, but more complexity, to the design process. In a given application, parameter studies are needed to find the optimal absorber parameters, which must include both performance and stability considerations. The methods demonstrated in this chapter provide the tools required for such investigations. In summary, the results of this dissertation have contributed to the design, analysis, and new applications of CPVA systems. While these results have not been experimentally investigated, previous studies have demonstrated the validity of models similar to those used herein, and therefore one has high confidence that the results are of practical utility. As a suggestion as future work, the techniques developed during this study can be combined in a simulation and design iteration environment. The HBM method proved to be quite powerful, and with proper design of the solver architecture can be extended to more complicated models for the “rotor”. This extention might include implemention of piston-crank dynamical models with periodic cylinder pressure excitation as well as simple drivetrain models. Then the guidelines to tuning of absorbers with flexible rotors can be utilized in this extended model structure. As the increased performance rocking absorbers are also understood in detail, this knowledge will provide the engineer more ways to achieve vibration reduction targets at given design constaints. 82 BIBLIOGRAPHY 83 BIBLIOGRAPHY [1] Tsio360 crankshaft. http://www.aircraft-specialties.com/tsio360-crankshaft-653139/. Accessed: 2017-03-05. [2] M A Acar and S W Shaw. Application of the harmonic balance method to centrifugal pendulum vibration absorbers. In Special Topics in Structural Dynamics, Volume 6, pages 243–252. Springer, 2016. [3] Torsional failures. http://www.drvconsulting.com/4-Broken-Crankshaft.jpg. Accessed: 2017-15-03. [4] E S Taylor. Eliminating crankshaft torsional vibration in radial aircraft engines. Technical report, SAE Technical Paper, 1936. [5] B Lecointe and G Monnier. Downsizing a gasoline engine using turbocharging with direct injection. Technical report, SAE Technical Paper, 2003. [6] R S Winberg. Internal crankshaft vibration damper, February 22 2000. US Patent 6,026,776. [7] S W Shaw and B K Geist. Tuning for performance and stability in systems of nearly tautochronic torsional vibration absorbers. Journal of vibration and acoustics, 132(4), 2010. [8] S W Shaw, P M Schmitz, and A G Haddow. Tautochronic vibration absorbers for rotating systems. Journal of computational and nonlinear dynamics, 1(4):283–293, 2006. [9] B Geist, C S Barron, W F Resh, M I Hossain, M P Patyi, and S M Mashkevich. Pendulum vibration absorber on a crankshaft, August 26 2014. US Patent 8,813,604. [10] W K Wilson. Pratical Solution of Torsional Vibration Problems Chapter XXX. Chapman & Hall, London, 1968. [11] J P Den Hartog. Tuned pendulums as torsional vibration absorbers. S. Timoshenko 60th Anniversary Volume, pages 17–26, 1938. [12] D E Newland. Nonlinear aspects of the performance of centrifugal pendulum vibration absorbers. Journal of Manufacturing Science and Engineering, 86(3):257–263, 1964. [13] J F Madden. Constant frequency bifilar vibration absorber, August 19 1980. US Patent 4,218,187. [14] H H Denman. Tautochronic bifilar pendulum torsion absorbers for reciprocating engines. Journal of Sound and Vibration, 159(2):251–277, 1992. 84 [15] A S Alsuwaiyan and S W Shaw. Steady-state responses in systems of nearly-identical torsional vibration absorbers. Transactions-American Society of Mechanical Engineers Journal of Vibration and Acoustics, 125(1):80–87, 2003. [16] A S Alsuwaiyan and S W Shaw. Performance and dynamic stability of general-path centrifugal pendulum vibration absorbers. Journal of Sound and Vibration, 252(5):791– 815, 2002. [17] B J Olson and S W Shaw. Vibration absorbers for a rotating flexible structure with cyclic symmetry: nonlinear path design. Nonlinear dynamics, 60(1-2):149–182, 2010. [18] B J Vidmar, S W Shaw, B F Feeny, and B K Geist. Nonlinear interactions in systems of multiple order centrifugal pendulum vibration absorbers. Journal of Vibration and Acoustics, 135(6):061012, 2013. [19] B J Vidmar, S W Shaw, B F Feeny, and B K Geist. Analysis and design of multiple order centrifugal pendulum vibration absorbers. In ASME 2012 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference, no. DETC2012-71074, Chicago, Illinois, USA, 2012. [20] R J Monroe, S W Shaw, B K Geist, and A H Haddow. Accounting for roller dynamics in the design of bifilar torsional vibration absorbers. Journal of Vibration and Acoustics, 133(6):061002, 2011. [21] V Cochelin and C Vergez. A high order purely frequency-based harmonic balance formulation for continuation of periodic solutions. Journal of sound and vibration, 324(1):243–262, 2009. [22] S Gozen, S W Shaw, and C Pierre. Resonance suppression in multi-degree-of-freedom rotating flexible structures using order-tuned absorbers. Journal of Vibration and Acoustics, 134(6):061016–061016–7, 2012. [23] J J Thomsen. Vibrations and stability: advanced theory, analysis, and tools. Springer Science & Business Media, 2013. [24] B F Feeny. Lecture notes in me 863: Nonlinear oscillations, Spring 2014. [25] Manlab - an interactive path-following and bifurcation http://manlab.lma.cnrs-mrs.fr/. Accessed: 2017-03-05. analysis software. [26] O Thomas, A Lazarus, and C Touzé. A harmonic-based method for computing the stability of periodic oscillations of non-linear structural systems. In ASME 2010 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference, pages 883–892. American Society of Mechanical Engineers, 2010. [27] J A Sanders, F Verhulst, and J A Murdock. Averaging methods in nonlinear dynamical systems, volume 59. Springer, 2007. 85 [28] A G Haddow and S W Shaw. Centrifugal pendulum vibration absorbers: an experimental and theoretical investigation. Nonlinear Dynamics, 34(3):293–307, 2003. [29] C P Chao, C T Lee, and S W Shaw. Non-unison dynamics of multiple centrifugal pendulum vibration absorbers. Journal of Sound and Vibration, 204(5):769–794, 1997. [30] R Cortright and D Farnsworth. Torsional vibration damper, April 27 2006. US Patent App. 11/412,520. [31] Masaru K, Takeshiro S, and Hideo O. An experimental study of a torsional/bending damper pulley for an engine crankshaft. In SAE Technical Paper. SAE International, 05 1989. [32] S Yu. Development of dual mode engine crank damper. Technical report, SAE Technical Paper, 2003. [33] C Huegel and C Dinger. Rotary vibration damper with centrifugal force pendulum, June 3 2014. US Patent 8,739,523. [34] C Shi, R G Parker, and S W Shaw. Tuning of centrifugal pendulum vibration absorbers for translational and rotational vibration reduction. Mechanism and Machine Theory, 66:56–65, 2013. [35] C Shi and R G Parker. Vibration modes and natural frequency veering in threedimensional, cyclically symmetric centrifugal pendulum vibration absorber systems. Journal of Vibration and Acoustics, 136(1):011014, 2014. [36] S W Shaw, M A Acar, B F Feeny, and B K Geist. Modal properties of rotating shafts with order-tuned absorbers. In Topics in Modal Analysis I - Proceedings of IMAC XXXII A Conference & Exposition On Structural Dynamics, pages 181–189, 2014. [37] M Häßler, S Lehmann, and A Rusch. Fliehkraftpendel und kupplungsscheibe mit demselben, 2013. DE Patent 112,011,104,426. [38] J Mayet and H Ulbrich. Tautochronic centrifugal pendulum vibration absorbers: general design and analysis. Journal of Sound and Vibration, 333(3):711–729, 2014. [39] T Krause, E Kremer, and P Movlazada. Theory and simulation of centrifugal pendulum absorber with trapezoidal suspension. In Proceedings of the 10th International Conference on Vibration Problems, September 5-8, 2011, Prague, Czech Republic, pages 322–327, 2011. [40] G Labahn and M Mutrie. Reduction of elliptic integrals to legendre normal form. University of Waterloo, Computer Science Department, 1997. 86