POSITION-DEPENDENT DYNAMICS AND STABILITY OF MACHINE TOOLS by Mohit Law B.E., University of Pune, India, 2003 MSc., Michigan Technological University, USA, 2008 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in The Faculty of Graduate Studies (Mechanical Engineering) THE UNIVERSITY OF BRITISH COLUMBIA (Vancouver) August 2013 ? Mohit Law, 2013 ii Abstract Machine tool?s productivity and ability to produce a component of the required quality is directly influenced by its dynamic stiffness at the tool center point. Lack of dynamic stiffness may lead to unstable regenerative chatter vibrations which are detrimental to the performance. The chatter vibrations are influenced by the changing structural dynamics of the machine as the tool moves along the tool path, resulting in position-varying machining stability of the system. Evaluation of these varying dynamics at the design stage is a complex process, often involving the use of large order finite element (FE) models. Complexity and computational costs associated with such FE models limit the analyses to one or two design concepts and at only a few discrete positions. To facilitate rapid exploration of several design alternatives and to evaluate and optimize each of their position-dependent dynamic behavior, a generalized bottom-up reduced model substructural synthesis approach is proposed in this thesis. An improved variant of the component mode synthesis method is developed and demonstrated to represent higher order dynamics of each of the machine tool components while reducing the computational cost. Reduced substructures with position-invariant response are synthesized at their contacting interfaces using novel adaptations of constraint formulations to yield position-dependent response. The generalized formulation is used to evaluate the position-dependent behavior of two separate machine tools: one with a serial kinematic configuration, and another with hybrid serial-parallel kinematics. The reduced machine model is verified against full order models and is also validated against measurements by including joint characteristics in the model. The effects of position and feed-direction-dependent compliances on machining stability are investigated by using a novel position and feed-direction-dependent-process-stability performance criterion that evaluates the productivity of machine tools in its entire work volume. Parameters limiting the target productivity levels are identified and modified; and, the complete dynamics are rapidly re-analyzed using the developed models. Optimal design modifications are shown to increase productivity by ~35%. The proposed methods in this thesis enable efficient simulation of structural dynamics, stability assessment as well as interactions of the CNC and cutting process with the machine tool structure in a virtual environment. iii Preface Parts of this thesis have been published in peer reviewed journals and conference proceedings. The publications reported work carried out during my Ph.D. research under the supervision of Dr. Y. Altintas and Dr. A. S. Phani. A version of Chapter 3 has been published as Law, M., Phani, A. S. and Altintas, Y., 2013, Position-Dependent Multibody Dynamic Modeling of Machine Tools Based on Improved Reduced Order Models, ASME Journal of Manufacturing Science and Engineering, Vol. 135 [1]. The reduced model dynamic substructuring strategy to model position-dependency was formulated by me with significant inputs from both supervisors. I wrote most of the manuscript. Parts of Chapter 4 have been published as Mehrpouya, M., Law, M., Park, C., Park, S. and Altintas, Y., 2013, Joint Dynamics Identification on a Vertical CNC Machine, 2nd International Conference on Virtual Machining Process Technology, Hamilton, Canada [2]. The section on joint identification in the Chapter is based on collaborations with M. Mehrpouya, C. Park, and Dr. S. Park from the Department of Mechanical and Manufacturing Engineering at The University of Calgary, Canada. I was responsible for virtual modeling, conducting experiments, and writing parts of the manuscript. A version of Chapter 5 has been published as Law, M., Altintas, Y. and Phani, A. S., 2013, Rapid evaluation and optimization of machine tools with position-dependent stability, International Journal of Machine Tools and Manufacture, Vol. 68 [3]. I modeled the position-dependent dynamics of the machine, conducted experiments and wrote most of the manuscript. Contents of Chapter 6 have been published as Law, M., Ihlenfeldt, S., Wabner, M., Altintas, Y. and Neugebauer, R., 2013, Position-dependent dynamics and stability of serial-parallel kinematic machines, CIRP Annals - Manufacturing Technology, Vol. 62 [4]. This work builds on our previous work [1,3] and was a result of collaboration with Dr. S. Ihlenfeldt, M. Wabner and Dr. R. Neugebauer from the Fraunhofer Institute of Machine Tools and Forming Technology IWU, Germany. I extended the position-dependent multibody dynamic machine model to model the dynamic behavior of a hybrid serial-parallel kinematic machine, conducted all the experiments, and wrote the manuscript. iv Table of Contents Abstract ..................................................................................................................................... ii Preface...................................................................................................................................... iii Table of Contents......................................................................................................................iv List of Tables ......................................................................................................................... viii List of Figures .......................................................................................................................... ix Nomenclature .......................................................................................................................... xv Acknowledgements ................................................................................................................ xxi Chapter 1 Introduction .............................................................................................................. 1 Chapter 2 Literature Review ................................................................................................... 11 2.1. Overview ........................................................................................................... 11 2.2. Modeling Machine Tool Dynamics .................................................................. 11 2.3. Modeling Position-dependency in Machine Tools ........................................... 12 2.3.1. Serial Kinematic Machine Tools ............................................................... 12 2.3.2. Parallel Kinematic Machine Tools ............................................................ 14 2.4. Model Order Reduction for Machine Tools...................................................... 15 2.5. Modeling and Identification of Machine Tool Joints ....................................... 19 2.6. Evaluation and Optimization of Machine Tool Performance ........................... 22 Chapter 3 Modeling Position-dependency in Machine Tools................................................. 26 3.1. Overview ........................................................................................................... 26 3.2. Machine Tool Component Modeling ................................................................ 27 3.2.1. Modeling Structural Substructures ............................................................ 27 3.2.2. Modeling Machine Tool Spindles ............................................................. 29 3.2.3. Modeling Feed Drive Units ....................................................................... 30 3.3. Substructure Model Reduction ......................................................................... 31 v 3.3.1. Improved Component Mode Synthesis Method ........................................ 33 3.4. Reduced Model Substructural Synthesis with Constraint Formulations .......... 36 3.4.1. Multipoint Constraint Equation Formulation ............................................ 38 3.5. Numerical Methods to Handle Constraint Equations ....................................... 39 3.5.1. Penalty Method .......................................................................................... 40 3.5.2. Lagrange Multiplier Method ..................................................................... 40 3.5.3. Substructural Synthesis.............................................................................. 40 3.6. Application ? Modeling Three Axis Milling Machine ..................................... 42 3.6.1. Verification of the New Improved Reduced Order Models ...................... 45 3.6.2. Structural Assembly for Different Positions of the Substructures ............ 48 3.7. Summary ........................................................................................................... 52 Chapter 4 Modeling and Identification of Machine Tool Joints ............................................. 54 4.1. Overview ........................................................................................................... 54 4.2. Modeling Joints between Structural Substructures ........................................... 55 4.2.1. Idealizing Bolted Connections ................................................................... 56 4.2.2. Idealizing Bearings as Springs .................................................................. 56 4.2.3. Idealizing Connections between Guide-ways and Guide-blocks, and between Ball-screw and Nut ............................................................................................ 57 4.3. Validating Model Predicted Response at Spindle Nose with Measured Response...............................................................................................................................57 4.3.1. Comparison of Model with Rigid and Spring Connections ...................... 58 4.3.2. Experimental Modal Analysis ................................................................... 59 4.3.3. Comparison of Measured and Updated Model Response ......................... 60 4.4. Finite Element Models of the Tool and Tool-holders ....................................... 62 4.5. Identifying Joints between Substructures ......................................................... 63 4.5.1. Inverse Receptance Coupling .................................................................... 63 vi 4.5.2. Joint Identification between the Tool and Tool-holder ............................. 66 4.5.3. Joint Identification between Spindle and the Tool-holder ......................... 71 4.6. Summary ........................................................................................................... 74 Chapter 5 Rapid Evaluation and Optimization of Serial Machine tools with Position-dependent Dynamics and Stability.......................................................................................... 76 5.1. Overview ........................................................................................................... 76 5.2. Generalized Substructural Formulation for Position-dependency .................... 77 5.2.1. Substructure Model Reduction .................................................................. 78 5.2.2. Substructural Synthesis.............................................................................. 79 5.2.3. Tool Point Response Predictions with Receptance Coupling ................... 80 5.3. Model Verification and Validation ................................................................... 81 5.3.1. Contrasting the Substructurally Synthesized Reduced Model with Full Model................................................................................................................................81 5.3.2. Comparing Model Predicted Response and Measured Response ............. 85 5.4. Position-dependent Dynamic Response at Tool Point ...................................... 87 5.5. Evaluation of Machine Tool Performance ........................................................ 90 5.5.1. Material Removal Rates and Position-dependent Stability ....................... 90 5.5.2. Oriented FRFs and Feed-direction-dependent Stability ............................ 90 5.5.3. Process Stability and Modes Limiting Productivity .................................. 93 5.6. Structural Design Modifications to Meet Target Productivity ......................... 95 5.6.1. Improved Machine Tool Performance ....................................................... 97 5.7. Summary ........................................................................................................... 99 Chapter 6 Position-dependent Dynamics and Stability Analyses of Serial-Parallel Kinematic Machines ............................................................................................................................... 101 6.1. Overview ......................................................................................................... 101 6.2. Position-dependent Multibody Dynamic Machine Model .............................. 102 vii 6.2.1. Substructure Model Reduction ................................................................ 104 6.2.2. Modeling Parallel Kinematics ................................................................. 104 6.2.2.1 Response Comparisons for Parallel Strut Modeled with Solid and Beam Elements.........................................................................................................................106 6.2.3. Reduced Model Substructural Synthesis ................................................. 108 6.3. Model Verification and Validation ................................................................. 109 6.4. Position-dependent Dynamics ........................................................................ 111 6.5. Evaluation of Position and Feed-direction-dependent Stability ..................... 113 6.6. Summary ......................................................................................................... 114 Chapter 7 Conclusions and Future Research Directions ....................................................... 116 7.1. Conclusions ..................................................................................................... 116 7.2. Future Research Directions ............................................................................... 118 Bibliography ......................................................................................................................... 120 viii List of Tables Table 3.1 Division of DOFs for the substructural components...............................................43 Table 3.2 Top 15 mode subsets for standard modal reduction (SMR) based method and mode indicator function (MIF) based sorting for the column substructure, where is the natural frequency..................................................................................................................................44 Table 3.3 Comparison of position-dependent behavior for the improved reduced model, where is the natural frequency, is the dynamic stiffness....................................................................................................................................51 Table 4.1 Modal Parameters obtained from measurements on FADAL 2216, where is the natural frequency, is the dynamic stiffness, and ? is the damping ratio..........................................................................................................................................60 Table 4.2 Comparison of dominant low frequency modes for model predicted response (with/ without joint effects) with measured response at the spindle nose, where is the natural frequency, is the dynamic stiffness, is the error in natural frequency prediction and is the error in dynamic stiffness prediction....................................................................61 Table 5.1. Division of DOFs for individual substructures.......................................................79 Table 5.2 Comparison of dominant low frequency modes for model predicted response (with/ without joint effects) with measured response at the TCP, where is the natural frequency, is the dynamic stiffness, is the error in natural frequency prediction and is the error in dynamic stiffness prediction.......................................................................................87 ix List of Figures Figure 1.1 Comparison of the traditional design process and the virtual design process........................................................................................................................................2 Figure 1.2 Flexible multi-body model of a machine tool..........................................................4 Figure 1.3 Typical computational effort required for position-dependent machine tool response analyses.......................................................................................................................5 Figure 1.4 Flow chart of thesis project ? Virtual integrated position-dependent process-machine interaction approach to design machine tools with targeted productivity...................8 Figure 2.1 Deviations in response due to unmodelled joints...................................................19 Figure 2.2 Process-machine interactions.................................................................................22 Figure 2.3 Block-diagram representation of chatter in milling with chip modulations...........23 Figure 3.1 Reduced model dynamic substructuring................................................................26 Figure 3.2 Substructural CAD and FE models for the column................................................28 Figure 3.3 Convergence checks to decide on the order of the substructures...........................29 Figure 3.4 Cross-section of Spindle unit and its FE model.....................................................30 Figure 3.5 Components constituting the feed drive model......................................................31 Figure 3.6 Column substructure represented only its interface DOFs.....................................32 x Figure 3.7 Flow chart for determining mode cut-off number..................................................34 Figure 3.8 Substructural assembly by enforcing continuity constraints, (a) compatible substructures, (b) incompatible substructures..........................................................................36 Figure 3.9 Overall flow chart showing a sequence of operations to obtain substructurally synthesized position-dependent response................................................................................41 Figure 3.10 Substructural assembly of the spindle-spindle housing substructure with the column substructure through constraint formulations.............................................................42 Figure 3.11 Normalized relative frequency difference (NRFD) and modal assurance criterion (MAC) comparison for standard component mode synthesis scheme and iterated improved component mode synthesis scheme ? for the column substructure.........................................46 Figure 3.12 Normalized relative frequency difference (NRFD) and modal assurance criterion (MAC) comparison for standard component mode synthesis scheme and iterated improved component mode synthesis scheme ? for the spindle-spindle housing substructure..............................................................................................................................47 Figure 3.13 Comparison of TCP FRFs for solution with different constraint formulations; assumed uniform damping of for all modes.............................................................48 Figure 3.14 Comparison of TCP FRFs for solution with different numerical methods; assumed uniform damping of for all modes.............................................................49 Figure 3.15 Comparison of full order model and reduced order model TCP FRFs at three different positions: top position (top), mid position (middle), and bottom position (bottom); assumed uniform damping of for all modes.............................................................50 xi Figure 4.1 Two-stage substructural synthesis of the machine tool. Tool (Substructure I) is synthesized with the tool-holder response (Substructure II); and, the synthesized tool-tool-holder combine is subsequently coupled to the response of the synthesized substructural machine model (Substructure III)............................................................................................55 Figure 4.2 Response comparisons at spindle nose for machine with rigid connections and spring connections in X and Y directions; assumed uniform damping of for all modes.......................................................................................................................................58 Figure 4.3 Experimental test setup for modal analysis on FADAL 2216 ? three axis vertical machining centre......................................................................................................................59 Figure 4.4 Measured response comparisons at spindle nose with model predicted response for rigid connections and spring connections in X (top) and Y (bottom) directions.....................60 Figure 4.5 Substructures coupled with a joint.........................................................................63 Figure 4.6 Free-free test setup for the tool-tool-holder combination.......................................67 Figure 4.7 Schematic of different tool and tool-holder assemblies.........................................67 Figure 4.8 Flow chart showing procedure for joint identification, and validation between the tool and the holder...................................................................................................................68 Figure 4.9 Identified translational joint FRF, httJ, between the tool and the tool-holder.......................................................................................................................................69 Figure 4.10 Identified rotational joint FRF, hrrJ, between the tool and the tool-holder.......................................................................................................................................69 Figure 4.11 Direct FRFs for the tool-tool-holder with 50 mm cylinder (tool), G11_50mm........70 xii Figure 4.12 Direct FRFs with tool-tool-holder with 90 mm tool, G11_90mm.............................70 Figure 4.13 Schematic of spindle, holder, tool assemblies......................................................71 Figure 4.14 Flow chart showing procedure for joint identification, and validation between the spindle nose and the tool-tool-holder combine........................................................................72 Figure 4.15 Identified translational Joint FRF between Spindle and tool-holder, httJ.............................................................................................................................................73 Figure 4.16 Direct FRF at TCP with 90 mm tool, G11_90mm.....................................................73 Figure 5.1 Flow chart of a virtual integrated position-dependent process-machine interaction approach for designing milling machines ensuring targeted productivity...............................76 Figure 5.2 Two stage substructural synthesis of the machine tool. Tool-tool-holder response (Substructure I) is coupled to the position-dependent response of the synthesized substructural machine model (Substructure II)........................................................................78 Figure 5.3 Full order and reduced order model TCP FRFs comparisons in X and Y directions..................................................................................................................................82 Figure 5.4 Normalized relative frequency difference (NRFD) and modal assurance criterion (MAC) comparison between the full order and reduced order machine model......................83 Figure 5.5 Experimental setup to measure TCP FRFs on machine.........................................85 Figure 5.6 Comparison of model predicted response (with/ without joint effects) with measured response...................................................................................................................86 xiii Figure 5.7 Comparison of reduced model TCP FRFs at three different tool positions: top, mid, and bottom.......................................................................................................................88 Figure 5.8 (a) Full order model (b) Reduced model with only interface DOFs, (c) Mode shape of column bending mode in YZ plane; (d) Mode shape of column bending mode in XZ plane.........................................................................................................................................89 Figure 5.9 Orientation of machine mode ( ? ) in the machine tool axes ( ) and the feed axes ................................................................................................................91 Figure 5.10 (a) Feed directional stability at three different tool positions: top, mid, and bottom, and, (b) Corresponding chatter frequencies as a function of feed direction...................................................................................................................................94 Figure 5.11 (a) Original, Modified, Optimized, and Smoothened Column models; (b) Comparison of X direction TCP FRFs for original and modified machine models; and, (c) Comparison of Y direction TCP FRFs for original and modified machine models......................................................................................................................................96 Figure 5.12 Feed-directional stability comparison for original and modified (optimized) machine model at two different tool positions: (a) bottom, (b) mid, (c) top...........................98 Figure 6.1. Hybrid serial-parallel scissor kinematic machine................................................102 Figure 6.2. Kinematic scheme for the scissor kinematics......................................................103 Figure 6.3. Oriented parallel strut modelled with Timoshenko beams..................................105 Figure 6.4 Solid Model of the parallel strut...........................................................................107 xiv Figure 6.5 Comparison of free-free response at free (connection) ends of the parallel strut modeled with different element types....................................................................................107 Figure 6.6 Overall flow chart showing a sequence of operations to obtain substructurally synthesized position-dependent response for the hybrid machine.........................................109 Figure 6.7 FRF comparisons between measured, reduced and full model............................110 Figure 6.8 Measured and predicted position-dependent behaviour of the machine in the X2Z plane.......................................................................................................................................112 Figure 6.9 Position-varying feed-direction-dependent stability over the whole work volume....................................................................................................................................114 xv Nomenclature Axial depth of cut Feed-direction-dependent critically stable axial depth of cut Coefficients of the characteristic equation A Cross-sectional area cA, cB Represent the locations on substructures A and B that are connected through the joint section xc Linear translational damping coefficient c? Rotational damping coefficient Displacement operator E Modulus of elasticity Error in natural frequency prediction Error in dynamic stiffness prediction f Total number of degrees of freedom (DOF) Natural frequency Natural frequency of the full order model Natural frequency of the reduced order model f Generalized force vector F(t) Cutting force as a function of time, t FiS Represents the vector of force and moment at location i on each substructure S FSJ Vector of force and moment in the joint section at the connecting locations to substructure S Cutting forces in the machine tool principal directions Cutting forces in the feed plane Represents the assembled structure?s frequency response function matrix h Represents the displacement-to-force receptance h(t) Dynamic chip thickness as function of time, t xvi httJ Translational frequency response function at the joint hrrJ Rotational frequency response function at the joint Static chip thickness H Frequency response function JH Frequency response function matrix of the joint Imaginary operator I Identity matrix k Interface node Bolt stiffness Dynamic stiffness xk Translational spring stiffness k? Rotational spring stiffness K Stiffness matrix Tangential cutting coefficient Radial cutting constant l Represents the displacement-to-couple receptance Number of modes of interest Bolt length m Number of constraint equations ? Structural mode M Moment Mass matrix n Substructure number nA, nB The internal locations on substructures A and B N Total number of frequency points Total number of elements Number of teeth o Represents the rotation-to-couple receptance p Total number of exterior DOFs being represented by the condensation node xvii P Number of component modes retained q Reduced order DOFs Generalized force/ couple vectors. Counter for number of modes Vector from the condensation node to the node corresponding to the interface node k R Subscript corresponding to reduced model matrices Real part Generalized receptance matrix that describes both translational and rotational component behavior Receptance matrix for the joint Rotation matrix Rotational operator to rotate element about the X axis Rotational operator to rotate element about the Z axis Represents the two-stage rotational operations s Number of frequency points S Substructure t Represents the rotations-to-force receptance T Tooth passing period Transformation matrix u Optimum response location Feed directional plane Displacement vector Condensation node displacement vector Displacement vector corresponding to exterior DOFs Displacement vector corresponding to interior DOFs Reduced displacement vector v Optimum driving point location Weight factor for DOF k y(t) Vibration marks at the inner modulation xviii y(t-T) Vibration marks at the outer modulation x, y, z Parameters pertaining to the x, y, and z directions (Cartesian) ,Sn cx Translational displacement vector for substructures A and B at locations nA, nB, cA, cB Generalized displacement/ rotation vectors Xi Displacement vector of the assembled structure Principal machine directions Penalty number Rotation angle about the X axis Penalty number matrix (diagonal) Orientation vector of the condensation node Matrix of directional factors Angles subtended by the platform Orientation of the mode ? with respect to the principal machine tool directions Rotation about the Z axis Lagrange multiplier Eigenvalues of the characteristic equation Real eigenvalues Imaginary eigenvalues Natural frequency [rad/sec] Chatter frequency Natural frequency corresponding to mode r [rad/sec] Mode shape vector for full order model Mode shape vector for reduced order model Eigenvector corresponding to the interior DOFs Transfer function matrix Oriented directional matrix 2D transfer function matrix at the tool point in the feed plane 2D transfer function matrix at the tool point in the principal plane xix Rotations about x, y, and z directions respectively Angular orientation ,Sn c? Rotational displacement vector for substructures A and B at locations nA, nB, cA, cB Damping ratio Instantaneous position of substructure Increments by which the substructures are adjusted to solve for the position-dependent response Acronyms CAD Computer aided design CMS Component mode synthesis CNC Computer numerical control DOF Degree of freedom FE Finite element FRF Frequency response function Improved CMS Iterated improved CMS IRC Inverse receptance coupling IRS Improved reduction system Iterated improved reduction system KE Kinetic energy MAC Modal assurance criterion MIF Mode indicator function MPC Multipoint constraint NRFD Normalized relative frequency difference ODP Optimum driving point RAM Random access memory RC Receptance coupling TCP Tool center point xx Mathematical Operators Includes, or is a superset of Is an element of Subset of Product over ? ? Norm of For each xxi Acknowledgements I would like to express my sincerest gratitude to my supervisor, Prof. Yusuf Altintas for his continuous encouragement and support all through my PhD years. His ideas and vision have helped guide, shape and inspire my research career. From him, I have come to better appreciate that the best science is that which is technically elegant and practically useful. I remain grateful to him for having given me this opportunity to be part of an extraordinarily talented and professional team at the Manufacturing Automation Laboratory. I am also deeply grateful to my co-supervisor, Dr. A S Phani for his invaluable support and mentorship. I owe much to him for making me appreciate, as he would put it ??that the devil is in the details?. From him, I have come to better understand that the best model of a cat is another, or preferably the same, cat. I would like to thank the NSERC CANRIMT Research Grant for the financial support provided for this research; and, for also making possible my research internship at the Fraunhofer IWU in Germany. I wish to also thank Dr. S Ihlenfeldt for his guidance and support during my stay in Germany; and the opportunities he made possible. I am also indebted to Prof. J Sutherland for his support in nudging me into applying to UBC to work under the supervision of Prof. Altintas. I remain grateful. My PhD years were shaped by interactions with many an exceptional people at the Manufacturing Automation Lab., and at UBC ? whom I have been fortunate to befriend and learn so much from. Thank you. All of my family ? thank you for having given me the space and support to wander in search. Medha ? my partner through all this and more, thank you. Many years ago, when I was younger still; I had read that many years later, my life would be in places that I cannot imagine still. Now, many years later, my life is in places that I would never have imagined possible. What has brought me here and will carry me further still are the people that I have met along the way. Thank you all. 1 Chapter 1 Introduction The rapidly evolving nature of the global manufacturing environment makes it necessary for machine tool manufacturers to respond to the dynamically changing market and product requirements in a timely manner. Moreover, as manufacturing paradigms continue to gravitate towards agility, rapid product output, shorter product development cycles and higher productivity requirements [5], machine tool builders can no longer afford traditional time-and-cost-intensive sequential product development cycles. Evaluation and optimization of machine tools is a complex process. Most concepts are evaluated against existing machines drawing on past design experiences, with little or no incorporation of formalized analysis methods within the design process. This is due to the fact that formal methods are too time-consuming to apply to more than one or two concepts. This results in an evaluation process where superior concepts are often overlooked or dismissed, resulting in inferior machine tool or costly modifications to a design late in the design cycle. Moreover, since up to 80% of the total life-cycle cost of a product is set during the conceptual stage of design [6], it becomes imperative that a systematic approach to design and development be followed so as to minimize the need for time consuming physical prototyping. ?Virtual prototyping? is a promising approach wherein the traditional time-and-cost-intensive sequential development process is replaced by an iterative process in which the machine is rapidly and iteratively analyzed and re-designed until performance requirements are met, see Figure 1.1. The increase in speed of analyses allows for an increased number of concepts to be evaluated and optimized for a given set of design (performance) constraints resulting in a higher quality concept while simultaneously reducing the whole product development time and cost significantly [7]. 2 Figure 1.1 Comparison of the traditional design process and the virtual design process [7] Performance of machine tools may be characterized by several different factors like machining speed, accuracy and reliability. However, the performance metric of main interest at the machine tool design stage is its achievable productivity characterized by material removal rates. This is governed by interactions between the cutting process and the structural dynamics of the machine. Lack of sufficient dynamic stiffness at the tool center point (TCP) may lead to unstable regenerative chatter vibrations. Chatter vibrations that are influenced by the machine tool?s dynamic stiffness are detrimental to the performance and integrity of entire machine tool system. These vibrations result in: poor surface quality, accelerated tool wear, damage of work piece and machine structural elements, and, ultimately limit productivity. Machine tool?s dynamic stiffness is the product of the modal stiffness and damping of a particular mode and is estimated as: ; where is the modal stiffness and is the modal damping. Machine tool stability is governed both by the tool-tool-holder and spindle dynamic characteristics and by the modal properties, i.e. the eigenfrequencies and mode shapes of the machine tool structural elements. These modal properties in turn are a function of the changing structural dynamics of the machine as the tool moves along the tool path in the machine work volume. The global lower frequency structural modes that belong to the 3 major machine tool components exhibit stronger position-dependent behavior as compared to the more application specific and local higher frequency spindle and tool-tool-holder modes. The position-varying dynamic response at the TCP results in position-varying stability of the system. Knowledge of the machine tool structure?s changing structural dynamics can be used to predict the cutting performance (i.e. chatter stability and geometrical accuracy), and to choose optimal values for the spindle speed and stable depths of cuts along the tool path. For large machine tools with large work volumes, the position-varying stability may require planning dynamically changing machining trajectories ? which pose its own set of challenges; or, alternatively and more conservatively, it may result in selection of cutting parameters below the lowest of all possible stability thresholds, thereby resulting in a slower material removal process. Furthermore, the position-varying structural vibrations between the TCP and servo drives also limit the positioning speed and accuracy of a machine tool. The objective of the machine tool designer is to maximize the dynamic stiffness between the TCP ? workpiece and the TCP ? drive motors while keeping the machine mass light for high speed positioning and high-productivity machining. Moreover, high-performance machine tools must also have consistently high dynamic stiffness over the entire work volume to mitigate position-dependent stability in order to minimize dynamically changing stable machining trajectories. Evaluation and optimization of a machine tool?s performance to deliver desired dynamic behavior over the entire work volume requires rapid assessment of several design alternatives in a virtual environment before eventual prototyping. Machine tool design and analyses involves the modeling of several subsystems by including: the machine tool behavioral model, the topological model, the kinematics model, control model, geometric error model, spindle dynamics, thermal error model, process model, workpiece attribute model and the fixture model [8]. A comprehensive treatment of all such models is beyond the scope of this present work ? which is mainly concerned with efficient position-dependent multibody dynamic modeling of machine tools, including characterizing the process-machine interactions via a stability model. Machine tools exhibit position-dependency because of varying boundary conditions as the tool moves along the tool path in the machine work volume. This motion may generally be achieved with a combination of linear axes motions for the case of machine tool with 4 serial kinematics or with change in the orientation, and/ or length of the parallel struts supporting the tool/ workpiece platform for the case of a machine tool with parallel kinematics. Alternatively, tool motion may also be achieved with a combination of serial and parallel motions for machine with hybrid kinematics. Modeling position-dependency for each of these machines presents its own set of unique challenges. This thesis develops generalized flexible multibody dynamic machine tool models to predict the position-dependent dynamic behavior for all configurations of machine tools. To facilitate such analyses, virtual machine tool models are treated as a collection of interconnected rigid and deformable bodies, as shown in Figure 1.2. Such a flexible multibody approach to machine tool modeling in which all the major machine elements are modeled with Finite Elements (FE) provides a comprehensive and complete description of machine tools, as both large movements of the bodies and structural deformations are modeled at the same time. Figure 1.2 Flexible multi-body model of a machine tool The FE method for machine tool design and analyses has proved useful for subsystem level design and optimization. However, response analyses of full machine FE models which are typically on the order of 1,000,000 degrees of freedom (DOF) or more is computationally costly and can take up significant portion of the total computational effort required for the design and analyses of machine tools [9]. Moreover, modeling position-dependency in such 5 large order FE models requires cumbersome adaptive/re-meshing strategies, making it time consuming in practice. Computations times up to several days is not uncommon to evaluate the changing structural dynamics over the complete work volume - as shown in Figure 1.3 for a typical machine tool with 1,000,000 DOFs. Figure 1.3 Typical computational effort required for position-dependent machine tool response analyses Position-dependency can alternatively be modeled using co-simulation techniques in which FE solvers are coupled to commercial multibody simulation codes [10]. However, these co-simulations are time prohibitive, and are limited to modeling machine tools as a combination of rigid-flexible bodies in contact; hence, are not completely able to capture the effect of component level flexibilities on the overall TCP response. To facilitate a computationally efficient modeling alternative to model position-dependency of machine tools in place of the presently used time consuming full FE models by the designers, a bottom-up approach based on substructural synthesis of flexible reduced models is proposed in this thesis. Main machine components are represented by their position-independent reduced substructural models which are subsequently synthesized to obtain position-dependent response. This bottom-up approach, also known as dynamic substructuring is a way to obtain the structural dynamic behavior of large and/or complex structures by dividing them into several smaller, simpler substructures (or components) for which the dynamic behavior is generally easier to determine. The dynamics of the total structure are then obtained by assembling the dynamic models of the individual components at a desired position. Since substructural reduction is carried out offline prior to synthesis, the 6 only serious computational effort is that required to solve the reduced synthesized equations of motion, thereby leading to significant computational savings. The substructural synthesis is akin to a generic component mode synthesis (CMS) approach [11], wherein each substructure modeled with FE is reduced to be represented by a set of interface (exterior) DOFs complemented by a set of generalized (modal) DOFs corresponding to the interior DOFs. For reduced models to have a manageable size, while retaining the essential dynamic characteristics of the full model, a major challenge with the CMS approach is to identify the number of modes from the individual components to be retained. To address this, a novel method to judiciously select only the significant modes that simultaneously represent higher order dynamics and keep the size of the reduced model to a minimum is proposed in this thesis. Furthermore, a generalized framework to incorporate mesh incompatibility at substructural interfaces is presented by describing displacement compatibility between the adjacent reduced substructures with sets of algebraic equations. These constraint equations are updated to account for relative motion; thereby simulating the position-dependent response. The accurate prediction of dynamic performance of machine tools at the design stage remains a challenge due to the difficulty in modeling joints within the FE environment. Joint characteristics depend on parameters like contact surface conditions, friction and damping; and since about 60% of the total dynamic stiffness and about 90% of the total damping in a machine tool structure originates at the joint [12]; joints must necessarily be considered in the model. Detailed and accurate modeling of the joint characteristics in machine tools is a separate and important research issue, and has been the subject of many past and on-going research investigations [13?15]. However, due to the variability in the multitude of mechanical interfaces in a typical machine tool, this becomes non-trivial. High fidelity joint modeling is not the main focus of this thesis; however, for the sake of completeness of the proposed multibody machine tool model, joint characteristics are obtained by using a two stage substructural synthesis approach. Joint characteristics between the tool, the tool-holder and the spindle are identified using a dynamic substructuring approach; and, other contacting interfaces are idealized as being connected by spring elements, damping for which is obtained by correlating model predicted response with measured response. 7 The current research aims to address the aforementioned challenges by developing a generalized and computationally efficient flexible multibody dynamic approach to modeling of machine tools to predict their position-dependent dynamic behavior. A generalized approach is sought which may be able to handle any kinematic configuration, including modeling machine tools with parallel kinematics. These virtual models will be further utilized to assess the position-varying machining stability of the system by modeling the process-machine interactions and generating productivity charts that characterize the productive cutting conditions within the entire machine tool work volume. The productivity charts quantify the speed-independent absolute stable depth of cut for every tool position within the machine work volume by also factoring in the dependence of machining stability on tool path/posture relative to the dynamically most compliant direction. These productivity charts in turn will form the basis for assessing machine tool performance by evaluating a given machine tool design concept for its ability to meet aggressive cutting conditions while machining/milling representative workpiece materials. The choice of workpiece material to be machined governs the selection of cutting/spindle speeds, which in turn influences the range of excitation frequencies. Modes around these excitation frequencies are typically prone to being self-excited during machining and may limit the productive cutting conditions. Milling of difficult-to-cut materials such as hardened steels is carried out at lower cutting speeds and tends to generate excitation frequencies below 200 Hz - the range corresponding to the global structural modes of the machine tool. Whereas, milling of free machining steels and light alloys generates excitation frequencies in the range of 200-1500 Hz, which correspond to the spindle and tool-tool-holder modes. Therefore, if the machine is to be conceived for milling hard materials, the stability of the process is dominated mainly by the modal parameters of the machine structure, whereas if the machine is to be conceived for machining light alloys, the stability is dominated by the modal parameters corresponding to the spindle and tool-tool-holder components [16]. Though the models developed in this thesis are generic, in that they may be used to assess the performance of machine tools for varied applications, the design evaluation and optimization will be limited to the case of machines envisaged for machining difficult-to-cut materials, and for the situation where the modes limiting the productivity correspond to the global structural modes of the machine which exhibit strong position-varying behavior. 8 Furthermore, the effects of varying workpiece dynamics, both in the case of large workpieces as well as thin-walled workpieces are ignored, and major position-varying behavior is assumed to be due to the machine alone. A flow chart of the overall scheme which facilitates rapid evaluation and optimization of a machine tool?s position-dependent performance is shown in Figure 1.4. For a given machine tool concept and defined aggressive cutting conditions the productivity charts are generated and the parameters limiting the target productivity levels are identified. These parameters are then subsequently modified (optimized) to meet design targets, using the virtual multibody dynamic multibody dynamic machine tool models that are developed. Figure 1.4 Flow chart of thesis project ? Virtual integrated position-dependent process-machine interaction approach to design machine tools with targeted productivity 9 Three design modules are proposed, namely: the improved reduction method, the generalized substructural synthesis algorithm, and the performance evaluation metrics through characterizing achievable material removal rates. Put together they form a platform for virtual simulation and optimization of machine tools. Though ?virtual prototyping? inherently assumes the absence of a physical prototype, to assess the quality of the virtual models developed in this thesis and to ensure that the models are physically relevant, the developed models are validated against measurements on similar available physical prototypes. Furthermore, since model validation necessarily requires the inclusion of joint characteristics in the model which involves measurement of damping from experiments, experimental modal analysis is carried out on the available physical prototypes, and this information is used in virtual model development. Since modeling/identification of machine tool joints requires damping to be measured, it is assumed that this information is available with the designer. Methods and models developed in this thesis are meant to aid the machine tool designer to efficiently explore several design alternatives during the initial stages of design and development. The developed models are judged adequate as long as they are reasonably able to capture behavior of the more detailed full order models. The developed models are meant to enhance engineering decision making prior to physical prototyping, and as such, the models will find limited applicability in establishing experimental guidelines based on model predictions. The rest of the thesis is structured as follows: first, a review of related literature is presented in Chapter 2, providing a backdrop based on the work of other researchers to evaluate methods presented in the subsequent chapters. Chapter 3 presents the generalized reduced model substructural synthesis approach to model position-dependency. At first, methods to model the substructural components are discussed, followed by dynamic substructuring. An improved substructural model order reduction scheme is presented and subsequently verified against full order models. The reduced substructural models are then synthesized at the desired configuration using novel adaptations of constraint formulations. In Chapter 4, a joint identification algorithm is presented to identify joints between the tool, the tool-holder and the spindle. Also, discussed are joint model updating techniques to 10 model the bolted connections as well as the rolling and sliding contacts. Model predicted response is validated against measured response for a typical three axis machine tool. The generalized model is extended to model and evaluate the position-dependent dynamic behavior of two separate machines with different kinematic configurations. First, a typical three axis milling machine tool is modeled and validated against similar available physical machine, which forms the contents of Chapter 5. Position and feed-direction-dependent-process-stability is evaluated to gauge the productivity of this machine in its entire work volume. Modes limiting productivity are modified and design optimizations are shown to increase productivity of the machine by as much as 35%. As a secondary application discussed in Chapter 6, position-dependency of a unique hybrid serial-parallel kinematic machine is modeled to evaluate its potentially superior dynamical behavior over traditional serial/ parallel kinematic machines. Since the dynamic stiffness inherently changes with position for machine tools with parallel kinematics, a rapid assessment of the influence of the change in dynamic stiffness on the machine capabilities is essential at the design and development stage for new machine tool architectures. This is greatly facilitated by extending the developed multibody dynamic machine tool model to evaluate this new concept by including a new method to model the simultaneously changing boundary conditions for the parallel struts undergoing translation and rotation. In the final chapter, concluding remarks and future research directions are discussed, followed by the bibliography. 11 Chapter 2 Literature Review 2.1. Overview The main objectives of this thesis are to model the position-dependent dynamics observed in machine tools and develop models that facilitate rapid evaluation and optimization of several design alternatives. This chapter reviews the literature relevant to the objectives of this thesis. At first, multibody dynamic models for machine tools and their extensions to model position-dependency are discussed in Sections 2.2 and 2.3. Since the dynamic substructuring approach proposed in this thesis relies on synthesizing reduced order models, a review of the relevant model order reduction schemes is presented in Section 2.4. Section 2.5 examines existing methods to joint modeling and identification. This is followed in Section 2.6 by a discussion on literature relevant to process-machine interactions which result in unstable machine tool chatter vibrations. Section 2.6 also includes a discussion on evaluation of machine tool performance based on their position and feed-direction dependent machining stability. 2.2. Modeling Machine Tool Dynamics Modeling of machine tool dynamics primarily involves obtaining the dynamic response at the tool center point (TCP). For this, the machine may be described by a simple lumped parameter model concept [17?19]; or the machine model may be derived from modal models [20?22]; or the machine may be modeled based on state-of-the-art finite element method representation [7,16,23?32]. Simplistic lumped parameter models [17?19] cannot completely capture the complexities of the entire machine tool system. On the other hand, modal models which are reconstructed from experimental modal analysis of physical machine tools [20?22] though accurate are restricted to characterizing only the measured machine; and may not easily be extended to model a different machine configuration during design process. The preferred contemporary approach is thus to describe the machine using finite element (FE) models. These models are efficient for subsystem level design analyses such as modeling of ball-screw feed-drive 12 systems [24?26], and spindles [27,28]. FE models for full machine analyses have also proven to be useful for structural modification based on process-machine interactions [12,16,23,28,29,31?33], as well as for control-structure interaction studies [30,33?36]. 2.3. Modeling Position-dependency in Machine Tools Position-dependency of machine tools is caused by position-varying boundary conditions as the tool moves along the tool path. Tool motion, made possible by different kinematic configurations, further influences the position-dependent behavior. This behavior, whether for machine tools with serial or parallel kinematics is generally modeled either by describing the machine by a FE model, or, by a flexible multibody dynamic model. Modeling position-dependency with large order FE models requires cumbersome adaptive/ re-meshing strategies for response analyses for each discrete position within the work volume, making it time consuming in practice. A multibody dynamic model on the other hand, traditionally uses a rigid body description of the machine components, and is appropriate to obtain quick and rough prediction of machine behavior; this, however does not capture the effect of body flexibilities on the overall response of the machine; nor does it completely account for the influence of kinematical configurations on position-dependency. 2.3.1. Serial Kinematic Machine Tools Recent methods to model position-dependency in machine tools with serial construction have focused on using co-simulation schemes in which FE solvers are coupled to standard flexible multibody dynamic analysis software [37?39]. A good overview of how co-simulation schemes are used for the modeling and simulation of the dynamic behavior of serial machine tools is provided by Reinhart and Wei?enberger in [10]. These co-simulation methods [10,37?39] have been shown to be effective for rigid-flexible body motion analyses and are less suited to the flexible bodies of a machine undergoing relative motion. Moreover, these co-simulations are time prohibitive, and are restricted to flexible bodies attached at geometrically fixed contact points, hence not suitable for flexible machine tool structural elements undergoing relative planar (translational) motion; which consist of multiple geometrically changing contact points that change as the tool moves from one position to another. 13 Other methods related to modeling the position-dependent dynamics for serial machine tools have largely centered on reducing the moving contact problem to that of a rigid body in contact with and moving over a flexible body with uniform nodal spacing, i.e. for bodies in contact that satisfy single set nodal compatibility conditions [25,40,41]. Position-dependency for multiple flexible bodies in contact was also modeled in [24,42,43] for a simplified single set nodal compatibility condition. Increasing modularity in the machine tool development process often requires substructures to be designed and modeled separately, resulting in different mesh resolutions at the contacting interface, that may not satisfy nodal compatibility conditions. Ensuring mesh compatibility during synthesis for such models which are simultaneously in contact over multiple nodes is difficult to guarantee, and if this condition was not necessary ? modeling time could be halved, as estimated in [9]. Substructures that are modeled independently and result in different mesh resolutions at the contacting interface during synthesis may be coupled by approximating the displacement field at the interface using a shape function formulation as proposed in [44?46]. These methods ensure displacement compatibility at the interface by using the shape function definitions of the adjacent elements being coupled, and can be quite cumbersome. To overcome these issues, Zatarain et. al. [9] proposed a dynamic substructuring approach to synthesize substructures with incompatible mesh types. A similar dynamic substructuring approach was also proposed by Fonseca in [34,47], however, it was limited to substructural synthesis for uniform nodal distributions at the coupling interfaces. An alternate and powerful frequency based approach to model position-dependency based on the so-called receptance coupling substructuring approach [48,49] in which the receptances of two independently modeled substructures are synthesized to obtain the combined response is also limited to the case of single set nodal compatibility and hence is not directly applicable to machine tools, in which substructures are simultaneously in contact over multiple nodes and may not have nodal compatibility at the interface. Different than the standard finite element methods, other mesh free (meshless) methods may satisfy the nodal compatibility criterion during substructural synthesis. However, these mesh free methods are reported to be slower than the standard finite element method [50], hence, are not suitable for the present application. Moreover, such mesh free methods are still 14 in their nascent stages of development and implementation, and as such, find limited applicability in the modeling of machine tools. 2.3.2. Parallel Kinematic Machine Tools Since high-performance machine tools are required to have consistently high dynamic stiffness over the entire work volume, nimbler parallel kinematic machine tools with higher dynamic stiffness capabilities are generally preferred. Since springs in parallel are known to be stiffer than springs in series; machine tools with parallel kinematics have a higher dynamic stiffness than machine tools with serial kinematics. However, as compared to machine tools with serial kinematics, machine tools with parallel kinematics exhibit stronger position-dependent dynamic behavior [51]. Though parallel kinematic machines have been studied extensively to identify their optimal kinematic configurations from a design point of view [52?55], their dynamics, which are known to vary considerably within the work volume, have been less investigated [54,56,57]. The co-simulation techniques that have been applied for serial machine tools [10,37?39] find limited applicability for machine tools with parallel kinematics; since these machines achieve tool motion by a change in the orientation, and/ or length of the parallel struts supporting the tool/ workpiece platform; resulting in flexible bodies that undergo simultaneous translation and rotation. To exploit the potentially superior dynamic behavior of machines with parallel kinematics, it is necessary to consider the position-dependency of the system at the design and development stage. A review of the relevant literature for machines with serial as well as parallel kinematics shows clearly that there exists a need for a generalized approach to model position dependency which treats large movements (translation and/ or rotation) of the bodies while including body flexibilities, i.e. structural deformations and being computationally efficient at the same time. To facilitate a computationally efficient modeling alternative to model position-dependency for both types of machine tools in place of the presently used time consuming full FE models by the designers, a bottom-up approach based on substructural synthesis of flexible reduced models is proposed in this thesis; including a novel method to model the parallel struts undergoing rotation. Also, a generalized framework tolerating mesh incompatibility at substructural interfaces is presented by describing displacement 15 compatibility between adjacent substructures with sets of algebraic equations, which are updated to account for relative motion; thereby simulating the position-dependent response. 2.4. Model Order Reduction for Machine Tools Modeling machine tool with finite elements results in very large order models, i.e. models with very high (~1,000,000) number of DOFs. These large order models are necessary to capture the topological and geometrical features in sufficient detail. However, response analysis for such large order models is computationally costly, making model order reduction a necessary consideration for efficient flexible multibody dynamic simulations. Reduction is typically characterized by a projection transformation from a higher order space to a lower order subspace with the goal of retaining original system characteristics with far lower storage and computational requirements. Model reduction methods to obtain the transformation matrix in structural dynamics can broadly be categorized as: physical methods, projection (hybrid) methods, generalized coordinate reduction methods, and modal methods. A detailed review on model order reduction techniques can be found in [15,58?60]; the following discussions treat only the most relevant methods to this thesis. Modal methods of model reduction, otherwise popularly also known as balanced reduction and/ or proper orthogonal decomposition methods [61?63] require converting the large scale second order ordinary differential equations obtained from the FE model to a linear first order state-space model. In the case of the balanced reduction method [63,64], the system is transformed to a basis where the states which are difficult to reach are simultaneously difficult to observe. Then, the reduced model is obtained simply by truncating the states which have this property. Though accurate, the modal methods are less suited to the present case for two reasons: they result in an increase in the size of the system from a second order to a first order state space decomposition; and because the reduced basis may not necessarily be associated with the physical coordinates of the interface DOFs, which is a requirement for the bottom-up substructural synthesis approach followed in this thesis. Since substructures are to be synthesized only at the contacting interfaces, it becomes necessary to reduce each individual substructure only to DOFs corresponding to the contacting interfaces. Hence, the reduction procedure involves eliminating a subset of DOFs from the displacement 16 vector, by partitioning it into the DOFs to be retained, i.e. the exterior /interface DOFs, which are in physical contact with the other substructure(s), and, the DOFs which are to be eliminated, i.e. the interior DOFs. Physical methods of model reduction [15,65?67] involve eliminating the interior DOFs such that the coordinates of the reduced model, i.e. the exterior DOFs represent only a subset of the full model. The simplest of this class, the static condensation method also known as Guyan reduction [65] is unable to represent higher order dynamics, since effect of the inertia terms corresponding to the interior DOFs are not considered. A modified form of the static condensation is the dynamic condensation method in which the inertia terms are considered at an appropriate initial frequency, the choice of which is non-trivial [15]. To overcome the shortcomings of the static and dynamic condensation methods, an improved reduction system (IRS) was proposed by O?Callahan [66]. IRS perturbs the static transformation by taking into account the inertia terms as pseudo-static forces. Though the IRS method offers a significant improvement over the static condensation method, the IRS still depends on the reduced mass and stiffness matrices obtained by static reduction. In order to minimize the error produced by this scheme, IRS could be extended to the iterated IRS method [67]. The iterated IRS algorithm was shown to converge to yield eigenvalues and eigenvectors of the full system. However, it results in a stiffer reduced system. Though the above discussed physical methods have their advantages, they are heavily dependent on the choice of the exterior DOFs, and yield problems when the exterior DOFs are clustered to one region of the substructure - as they are in the case of the exterior DOFs representing only the interface DOFs in the present case. The generalized coordinate reduction methods [59,68] offer an alternate form of model reduction which is independent of the selection of the exterior DOFs. Of these methods, the System Equivalent Expansion Reduction Process (SEREP) [68] results in ?exact? representation of the full order model up to a predefined frequency range. Though ?exact?, SEREP is computationally inefficient since it requires computation of the eigenvectors of the full order model system of equations. Alternatively, the Ritz vector method which is computationally more efficient that computing the ?exact? eigenvectors may be utilized. However, the dynamic equations of motion of the reduced model obtained from Ritz vector methods are generally coupled [59]. 17 Projection based model order reduction methods include methods such as the Component Mode Synthesis (CMS) and its many variants [11,60]. Most applications of CMS employ one of the two approaches: the fixed-interface-mode methods, or, the free-interface-mode-methods. The former, also known as the Craig-Bampton reduction method, is the more widely used form of the CMS methods. It employs fixed-interface component normal modes and constraint modes by retaining all physical boundary DOFs as independent generalized coordinates, greatly facilitating substructuring and hence forms the basis for the model order reduction scheme employed in this thesis. The constraint modes are similar to the static condensation method (Guyan reduction), ignoring any inertial contributions. The IRS based reduction was combined with by Koutsovasilis [69] to include inertia terms. Further discussions on definitions of constraint modes and component normal modes are addressed in Chapter 3. The effectiveness of the Craig-Bampton method is still heavily dependent on the definition of exterior DOFs and the number of component normal modes retained. Since the number of exterior DOFs is fixed by the definition of DOFs belonging to the contacting interface, we shift our attention to determining the number of component normal modes representing the interior DOFs which need to be retained. A major challenge with the CMS approach is to identify the number of mode sets from the individual components to be retained in order to capture the essential dynamics without increasing the model size. Usually, either the first few low frequency modes; or modes with an eigenfrequency up to 1.5-2 times the maximum excitation frequency are retained [70]; this however may increase the size of the reduced model. Moreover, it may not be sufficient to represent the flexibilities of substructures using only a single mode set for each component due to their position-varying boundary conditions [41]. An accurate representation of the full frequency spectrum while retaining computational efficiency requires answering the following question: how many and which modes of the subsystem need to be retained? Modal methods of mode selection [71,72] that are based on ?dominance measures? involve sorting and ranking eigenvalues of the system in order of their importance (dominance) also require decomposition into a linear state-space model, resulting in even larger system matrices, hence are not preferred here. Other methods [70,73,74] sort and retain the most significant modes based on the strength of interaction between the subsystem 18 and the main system. In the present situation, as one machine tool substructure moves over another, these methods [70,73,74] may result in identification of different significant mode sets for different positions, due the nature of the changing interactions. This is contrary to the objective of synthesizing substructures that have position-invariant response. A mode selection criterion was also proposed by Park in [75]. While being reliably able to select important modes, their top-down approach was based on partitioning of already meshed continuums into substructures, thus resulting in a conforming mesh at the interfaces. Hence, is different than the generalized bottom-up approach followed in this work such as to synthesize substructures which may have been modeled separately and may have non-conforming mesh distribution at their interfaces. Alternatively, suitable subspace dimension determination may be carried out using an adaptive refinement method as pursued in [76], wherein modes were included/excluded iteratively from each subspace depending upon how posteriori error estimates of each CMS subspace influenced the error in the synthesized reduced model. Though automated, it is computationally expensive and inefficient. To address the above mentioned issues, and to develop substructural reduced models that have position-invariant response that are able to represent the essential dynamics of the full system, it is essential to retain only the significant component modes of each substructure. To acheive this, a novel mode indicator function (MIF) based identification criterion is proposed in this thesis which identifies modes simulated from a frequency response function (FRF) constructed between the most controllable and observable locations within the subset of the interior DOFs. FRFs are preferred in this investigation because they provide complete information about a system?s dynamic behavior by including information about the modal stiffness, damping and natural frequencies of the system. The mode sets thus identified with the MIF from the simulated FRF are treated as significant. These modes are expected to represent higher order dynamics of the substructure while keeping the order of the reduced model to a minimum by spanning a much wider frequency range with fewer modes than would be required with standard CMS methods, as discussed in Chapter 3. Each reduced substructure is combined with the other substructures to obtain the combined response. During synthesis of the substructures, joints between the substructures also need to be modeled. 19 2.5. Modeling and Identification of Machine Tool Joints A typical machine tool includes many types of joints, each with different characteristics that affect the overall machine dynamic response differently. The bolted connections between structural members; the connections between the guide-block and rail and between ball-screw and nut; and the bearing supports ? all have varying degrees of influence on the low frequency TCP response; whereas the high frequency TCP behavior is more influenced by the interface connections between the tool and tool-holder, and between the tool-holder and the spindle respectively. Joint characteristics for each of these interface types depend on a multitude of parameters like preloads, contact surface conditions, bearing types, friction, and damping; information about which is seldom available at the design stage, making joint modeling non-trivial in the FE environment. The joints, if unmodeled, often result in deviations in the response characteristics between the virtual model and their corresponding physical prototypes, Figure 2.1. Figure 2.1 Deviations in response due to unmodelled joints Since the fundamental barrier to accurate predictive structural dynamical models is the variability of the mechanical interfaces, several accurate but complex models have been developed to approximate the bolted connections, and the connections between the guide-20 block and rail and between ball-screw and nut [32,77?85]. The contact stiffness at the contacting interfaces between the guide-block and the guide-rail, and between the ball-screw and nut depends on: the preload, the contacting surface preparation, the number of contact points between the ball grove and the raceway as well as the profile of the ball grove; and has been treated in sufficient detail in [32,77?80]. Nonlinearities in these translational joints due to preload and the nature of contact were also treated in [80?83] based on a parametric model and Hertzian contact mechanics. An interesting approach to joint modeling for these interface types was also proposed by Guo et. al. in [84], in which a FE model of the whole machine tool was used in conjunction with a lumped parameter machine model to obtain joint stiffnesses from dedicated detailed modal testing on the assembled machine tool structure. A detailed comparison of several methods to model bolted connections in machine tools with FE modeling was also discussed in [85]. Overall, though significant advances have been made in modeling these joint types, i.e. the bolted connections; connections between the guide-block and rail and between ball-screw and nut; and bearings, such high fidelity joint models necessarily require several dedicated experiments for validation accompanied by detailed FE modeling of the joint interfaces, making it difficult to extend for all such joints in a complete machine tool. Moreover, regardless of the joint type, damping must always be measured from experiments and in turn be used to update the FE model which is typically modeled without damping. Several studies have therefore focused to address identification of joint dynamics by updating the finite element model from experiments based on direct and/or indirect methods. The following discussions outline some of the most relevant methods to this thesis, and more comprehensive reviews may be found in [13,15,86]. Direct methods of finite element model updating [87,88] are based on updating the global assembled FE matrices in a single step resulting in new system matrices which may have little or no relevance to the physical test structure. The direct methods were superseded with sensitivity based methods [89] which required determination of the sensitivity of a set of updating parameters (typically the springs coupling different substructures) to differences in dynamic behavior between numerical/ analytical and experimental data. These methods which are iterative in nature are based either on measured modal data or measured FRFs, and require direct measurements at the joints and are thus impractical due to the inaccessibility of 21 the joints in assembled machine tool structures. Furthermore, all updating techniques pose difficulties in matching the order of the measured and simulated models; and also that measurement data is sometimes incomplete and that it may be corrupted by noise. An alternate to the direct joint identification technique is the force-state mapping technique [90,91]. In this technique, the force transmitted by the joint is represented as a function of the displacement and velocity of the joint, i.e. as a function of the full mechanical state of the joint. This allows rate-dependent effects associated with the joint to be identified. Though powerful, the force-state mapping technique necessitates several dedicated experiments, making identification very cumbersome. Hence, indirect methods of identification through the receptance-based approaches find favour. The receptance coupling (RC) approach couples experimentally or analytically obtained FRFs which are subsequently used to derive the response of the assembled structure. The RC method has been successfully employed to couple the receptances of different machine components and obtain FRFs at the TCP [92?95]. Effective deployment of the RC method requires measurements of the rotational FRFs, which are difficult to measure. Park and Altintas [93] used the RC method to extract the rotational FRFs of the tool and the tool-holder by performing two sets of measurements on the blank tool. Within the same family of the RC methods, the so-called ?inverse receptance coupling? (IRC) method has been employed to obtain the joint properties between different substructures by using the dynamics of the assembled structure along with the dynamics of the subcomponents [92,96]. In the IRC method, the dynamics of the assembled structure and subcomponents are employed to obtain the joint properties between different substructures. Since the IRC method requires measurements of the substructural dynamics in its un-assembled configurations, it is difficult to extend to the case of the entire machine tool systems. A machine tool is seldom available in its unassembled state to perform free-free experimental modal analysis on its substructures; and even if it is, measurements require time and cost prohibitive experimental platforms to be built. To incorporate joint characteristics that possess physical fidelity and that are sufficiently simple of form, such that they can be easily integrated into virtual machine tool models, a two-stage substructural synthesis approach is followed in this thesis. At first, response at the spindle nose is obtained from the substructurally synthesized reduced machine model by 22 idealizing joints between the various structural components of the machine as connected with linear springs; damping, for which is incorporated by correlating model predicted response with measured response. As a second step, response at the spindle nose is synthesized with the tool and tool-holder response by formulating an IRC approach which identifies the joint characteristics between these interfaces. The IRC method overcomes the barrier of identifying joint characteristics even where they cannot directly be measured. This, along with using the virtual machine tool model reduces the need to perform several dedicated measurements that were required in earlier studies to obtain translational, and, the more difficult to measure rotational frequency response functions. 2.6. Evaluation and Optimization of Machine Tool Performance The substructurally synthesized reduced machine models that are developed in this thesis allow for efficient prediction of the position-dependent dynamic response at the TCP. This further facilitates rapid investigation of the performance of the machine tool based on evaluating its position-varying chatter stability within the work volume. Chatter, an artifact of process-machine interactions results from the regeneration of chip thickness during the cutting process [97]. As shown in Figure 2.2, it is possible that during cutting, the cutting force signal has enough energy (frequency) content within it, that it excites one of the structural modes of the machine, resulting in modulation of the chip thickness. Figure 2.2 Process-machine interactions 23 A block diagram describing the system?s dynamics in the Laplace domain is shown in Figure 2.3, which includes a simplified 2 DOF representation of the milling system and a representation of chip modulations. The wavy surface being left on the workpiece by the cutting tooth in contact, known as the inner modulation, decreases the dynamic chip thickness, shown as h(t) in Figure 2.3, while the vibration marks left from the previous pass, y(t ? T) known as the outer modulation, increases the dynamic chip thickness; where T is the time period between two successive teeth. Therefore, when there is a phase shift between y(t) and y(t ? T), the dynamic chip thickness varies at the frequency of vibration and creates a vibrating cutting force F(t) which could amplify the vibration of the tool. If the energy diverted from the machining process by chip thickness variations is larger than the vibration damping capacity of the structure, the amplitude of the vibrations will grow until the vibration amplitude is large enough to make the entire dynamic system unstable, possibly resulting in tool breakage and damage to the workpiece and machine. Figure 2.3 Block-diagram representation of chatter in milling with chip modulations 24 Modeling these process-machine interactions for the purpose of predicting the onset of chatter in machining has been the subject of investigation for several decades. A good review of chatter may be found in [98]. The prediction of chatter stability requires the modeling of dynamic chip thickness, a cutting force model, workpiece material properties, cutting tool geometry, and most importantly, a dynamic model of the machine ? either measured or model based. Of these parameters, of paramount importance to the machine tool designer is accurate prediction of a machine?s tool point dynamic stiffness. Design optimization of structures to evaluate and subsequently enhance dynamic stiffness has been the focus of many past studies [17?20]. Early work by Tlusty [20] was based partially on a process-machine interaction approach to machine design in which chatter was proposed as a metric to evaluate the productivity of a typical turning machine. Reddy and Rao [17], and, Rao and Gandhi [18] proposed design optimization methods for minimizing the total weight of the structure under the constraints such as static rigidity, natural frequencies, and regenerative chatter stability. Yoshimura et. al. [19] used simplified structural models in order to minimize manufacturing cost of machine under constraints of machining accuracy, machining productivity, and local deformations of structural members. Rahman et. al. [22] also used stability charts to evaluate the effect of different materials for the machine base on the performance (productivity) of machine tools. These methods [17?20,22] though effective rely upon simplistic lumped parameter models or difficult to obtain modal models. More recently, Catania et. al. [21] proposed a hybrid model in which an experimental modal model of the machine was combined with analytical beam models for the tool to predict and evaluate machine performance via chatter stability. Recent studies [16,31,32,99,100] have also used more detailed finite element models of the machine to assess machine tool performance via a stability model. Since evaluation and optimization of machine tools is a complex process, there exists a need for a rapid assessment scheme of several design alternatives before eventual physical prototyping. Moreover, since the machine exhibits position-dependent dynamic behavior which results in position-dependent stability, there also exists a need to evaluate the machine performance within its active working volume. The computationally efficient substructurally synthesized reduced machine model developed in this thesis facilitates such evaluation. 25 It is well established that chatter stability depends on tool geometry, cutting conditions and on the machine?s dynamic response; whereas it is less known and consequently less investigated that chatter stability also depends on tool path/ posture relative to the dynamically most compliant direction. A multi-degree-of-freedom machine tool may chatter in any of its dominant modes. The effect of each mode is determined by its dynamic characteristics; and, whether or not the mode is aligned with the principal machine tool directions. Stability is further governed by the projections of these modes into the feed directions, which must also be factored in during evaluations of the attainable productivity within the entire machine work volume. Polar plots, which are a convenient measure to quantify the feed-direction dependent stability, were used in early studies by Tlusty in [20]. These plots chart the speed-independent absolute stable depth of cut as a function of feed orientation (0-360?). Tlusty [20] utilized these charts to determine the least compliant direction for different mounting orientations for the tool post in the case of a turning machine; which was later re-designed to increase its dynamic stiffness in the direction of the most complaint mode. Shamoto et. al. [101,102] proposed a somewhat similar stability index to optimize tool path/ posture to avoid chatter vibration in the case of ball-end milling with different tool inclinations and also in the case of turning to determine optimum tool feed angles ? much like in [20]. Zuliaka et. al. [16] also used the feed-direction dependent stability criterion to investigate effects of structural modification on feed-directional compliances in the case of a milling machine. Kilic and Altintas [103] also investigated the effect of unmatched spindle dynamics on feed-directional stability and suggested having matched spindle dynamics to have near uniform feed-direction dependent stability. In this thesis, the work presented in [103] is advanced to propose a generalized feed-direction-dependent stability criterion to characterize machine performance by its position and feed-direction dependent stability. This criterion also identifies parameters and components which limit target productivity levels. These components are subsequently modified to result in an improved machine concept that delivers target performance requirements. 26 Chapter 3 Modeling Position-dependency in Machine Tools1 3.1. Overview A new generalized bottom-up method is proposed in this chapter in which substructural elements of the machine tool that have position independent response are reduced independently and combined subsequently to obtain position-dependent dynamic response. A schematic overview of the proposed reduced model substructural synthesis formulation and its development in this chapter is shown in Figure 3.1. Figure 3.1 Reduced model dynamic substructuring A three axis vertical milling machine with serial configuration is selected as the representative machine to be modeled and analyzed. At first, from a given 3D Computer 1 A version of this Chapter has been published as Law, M., Phani, A. S. and Altintas, Y., 2013, Position-Dependent Multibody Dynamic Modeling of Machine Tools Based on Improved Reduced Order Models, ASME Journal of Manufacturing Science and Engineering, Vol. 135(2) [1] 27 Aided Design (CAD) model of the machine under investigation ? a representative three-axis vertical machining center; the machine is partitioned into its various substructural components. Following partitioning, each of the main substructural components of the machine is modeled in FE in Section 3.2. Substructural level model reduction will be carried out using an improved variant of the component mode synthesis method in Section 3.3. Reduced substructures are synthesized at the contacting interfaces by ensuring displacement compatibility by sets of algebraic constraint equations, which is treated in Section 3.4. Synthesis is carried out with two constraint formulations; followed by two numerical methods to handle the constrained equations of motion in Section 3.5. The chapter is concluded by successfully and efficiently simulating position-dependent TCP response for several different tool positions; and the model is also verified against full order model results. 3.2. Machine Tool Component Modeling Each of the main substructural components of the machine, i.e. the column, bed, table, cross-slide, spindle housing, spindle, and the three separate feed drives are all modeled with finite elements using ANSYS? [104]. After necessary convergence tests on FE models, the substructural system matrices are exported into the MATLAB? environment for further model reduction investigations. 3.2.1. Modeling Structural Substructures FE models for structural substructures have been generated from their respective detailed CAD models using 10-noded solid tetrahedral elements. These elements are preferred over other 8-noded solid brick elements, since unlike brick elements they do not cause any artificial shear locking [105]. The structural components, made of a grade of cast iron were assigned material properties as follows: modulus of Elasticity of 89 GPa; density of 7250 kg/m3; and, Poisson?s ratio of 0.25. The CAD model along with the FE model for one such structural component, the column, is shown in Figure 3.2. 28 Figure 3.2 Substructural CAD and FE models for the column Convergence tests were carried out to determine the order of the substructural finite element models. The h-method of refinement in which the size of the mesh is varied until results converge is employed in this work. This method, due to its ease in implementation is preferred over the p-method of refinement in which the order of the polynomial used for the shape function is changed while keeping the mesh size the same. Modal analyses investigations were carried out for each of the substructures with several different mesh sizes and the convergence rate for the first 100 natural frequencies were compared. The working definition of convergence criterion employed in this work is an error of up to 5% in in natural frequencies over the frequency range of interest for a given model size when compared with the natural frequencies obtained from the model having the highest number of DOFs. Figure 3.3 shows convergence checks for the first ten non-rigid body modes, for the example of the column substructure. For a 5% error in natural frequencies, a model size of 17103 DOFs is found adequate. If however, the convergence criterion was set as 2%, the model size required would be considerably higher at 102681 DOFs with significantly higher demand on computing resources. Mesh refinement was carried out using the ?auto? meshing option within ANSYS? which results in refinement of mesh only at locations with abrupt changes in cross-sections and geometry [104]. Hence, unlike structures with uniform mesh distributions for which the 29 lower frequency modes converge quicker than the higher frequency modes and the convergence is generally monotonic [106]; in the present case, the non-monotonic convergence of modes is thought to be due to the non-uniform mesh refinement procedure. Similar convergence checks to determine the size of the model for all other major substructures were carried out. Figure 3.3 Convergence checks to decide on the order of the substructures 3.2.2. Modeling Machine Tool Spindles The spindle assembly, shown in Figure 3.4, includes the tool-tool-holder, spindle shaft, spindle cartridge, bearings, spacers, drive pulley, and other accessories such as nuts and rotary couplings. This spindle assembly is modeled as described in detail in [27,107]. The spindle shaft, spindle cartridge, and the tool-tool-holder combination are all modeled with Timoshenko beam elements. Bearings are modeled as radial-axial springs, and other accessories are modeled as lumped mass elements. The tool-tool-holder-spindle interface connections are assumed to be rigid at this stage of the investigation; Chapter 4 deals 30 separately with modeling and identification of these joints. This comprehensive spindle assembly model has been previously validated against measurements in [27,107] and is integrated as a separate substructure coupled rigidly to the spindle housing. Figure 3.4 Cross-section of Spindle unit (top) and its FE model (bottom) [27] 3.2.3. Modeling Feed Drive Units Feed drive units for each of the three linear axes are modeled as separate subsystems. The mechanical elements constituting the feed drive model are shown in Figure 3.5. The ball-screws are modeled with Timoshenko beam elements. The translating unit (table, cross-slide or the spindle housing) are structural components modeled as described earlier in Section 3.2.1. Support bearings and the connection between the ball-screw and the nut are modeled as radial-axial springs [25,40]. The motor and other accessories are modeled as lumped mass elements. 31 Figure 3.5 Components constituting the Feed drive model [25] 3.3. Substructure Model Reduction Each of the main substructural machine components is reduced independently as discussed here. For any substructure under consideration, the undamped equations of motion are represented as: ? (3.1) where { } are the mass and stiffness matrices in real (physical) space; and is the force vector. In many cases, ; where represents the total number of degrees of freedom. The goal of model order reduction in structural dynamics is to find a low order subspace , to approximate the displacement vector in Eq. (3.1) such that: (3.2) where is the reduced displacement vector; and is the transformation matrix used for reduction; making the reduced structural matrices: (3.3) Since the reduction procedure involves eliminating a subset of DOFs from the displacement vector of the full model, , it is partitioned into the DOFs to be retained, i.e. the exterior /interface (E) DOFs, , that are in physical contact with the other substructures(s); and, the DOFs which are to be eliminated, i.e. the interior (I) DOFs, , as: [ ] { ? ? } [ ] { } { }. (3.4) 32 As an example of the interface DOFs to be retained, the exterior DOFs for the column substructure representing the guide-way top surfaces and the surfaces in contact with the machine base are shown in Figure 3.6. Figure 3.6 Column substructure represented only its interface DOFs Since the substructural assembly formulation is equivalent to a generic component mode synthesis (CMS) [11,46]; the standard form of the transformation matrix using the CMS - Craig and Bampton method, which contains the entire set of constraint modes corresponding to the exterior DOFs and a set of component normal modes corresponding to a subset of the interior DOFs is given as [11]: [ ] (3.5) where is a unit matrix; and, is the equivalent static (Guyan) transformation [65]. in Eq. (3.5) is obtained by traditionally retaining the first few P modes of the mode shape vector ; which is obtained by solving the eigenvalue problem corresponding to the interior DOFs of the form: [ ] (3.6) Though the standard form of the CMS transformation matrix in Eq. (3.5) reasonably approximates response, the quality and accuracy of the reduced model is limited by the quasi-static nature of the transformation as well as its dependence on the number of 33 component modes included in the transformation. Since the bottom-up formulation followed in this thesis is dependent on the substructures having position-independent response, it becomes necessary for the reduced substructures to be able to represent essential full model characteristics while simultaneously keeping the order of the reduced model to a minimum. 3.3.1. Improved Component Mode Synthesis Method To overcome the part-static nature of the transformation, an improvement as suggested by Koustaviils [69] is employed which includes the inertial terms as well; and, a novel method to judiciously select the most important component modes is proposed. 3.3.1.1 Addressing the quasi-static nature of the CMS transformation Advantageous characteristics of another model order reduction scheme, the so-called iterated ?improved reduction system? (IRSi) method based on [66,67,69] are included to overcome the quasi-static approximation. The new iterated improved component mode synthesis transformation matrix becomes [69]: [ ] (3.7) wherein the iterated IRS terms are [69]: ( ) (3.8) with [ ] (3.9) and (3.10) wherein . (3.11) The reduced matrices in Eq. (3.8-3.10), i.e. ; are obtained with appropriate projections onto reduced subspaces as in Eq. (3.3). 34 3.3.1.2 Component Mode Selection Criterion To identify how many and which of the P modes within the mode shape vector , i.e. within Eq. (3.7) are to be retained, the procedure as shown in the flow chart in Figure 3.7 is proposed. Figure 3.7 Flow chart for determining mode cut-off number Much like in Experimental Modal Analysis, where actuators are placed at locations to excite all the modes in the frequency range of interest, and, sensors are placed at locations where these modes achieve a high degree of observability; a similar methodology is used in this work, wherein the significant modes are treated as those identified with a mode indicator function from the simulated FRF constructed between point(s) of high controllability and observability within the subset of the interior DOFs. To determine the optimal excitation location, the ?optimum driving point(s) (ODP)? method [108] is employed since it offers a very ?natural? selection criterion based on the computed eigenvectors ( ). Point(s) with large mode shape amplitudes for the modes of interest are selected while avoiding the nodal points or any points near the nodal line. The ODP parameter for a DOF, v, is defined as [108]: 35 ?? ? (3.12) where is the number of modes of interest. The DOF that has the largest ODP parameter is considered as the possible excitation location. To determine the response location, the point with the maximum kinetic energy is chosen. The energy distribution method is employed here, as it involves computing the kinetic energy, which is a direct measure of the structural properties. The kinetic energy (KE) at the DOF for modes of interest is [109]: ? (3.13) where is the undamped eigenfrequency obtained by solving the eigenvalue problem in Eq. (3.6); and are the inertial terms corresponding to the interior DOFs obtained from Eq. (3.4). A FRF ( ) is simulated up to the frequency range of interest as: ? (3.14) where is the eigenvector at the optimal response location DOF, ; is the eigenvector component at the optimal driving point DOF, ; is the damping ratio; is the natural frequency; and, is the imaginary operator. To identify the number of modes in this FRF, which are to be treated as significant, a Mode Indicator Function (MIF) is employed which is defined at each frequency point, , for a total of N points as [108]: ? | || | ? ? ? (3.15) The transformation of Eq. (3.7) complemented with proper sorting and selection of the component modes to be retained, identified from Eq. (3.15) results in reduced order substructural models which are able to represent higher order dynamics of the substructure while keeping the order of the reduced model to a minimum by spanning a much wider frequency range with fewer modes than would be required with standard CMS methods. 36 3.4. Reduced Model Substructural Synthesis with Constraint Formulations Consider two substructures and as shown in Figure 3.8 which may represent any of the machine components with relative motion between them. These substructures are already reduced to their interface DOFs using the reduction scheme discussed in Section 3.3. Figure 3.8 Substructural assembly by enforcing continuity constraints, (a) compatible substructures, (b) incompatible substructures If the substructures are obtained through a partitioning of the FE mesh of the global structure, i.e. as in a classic top-down dynamic substructuring approach, they would be automatically conforming, as in Figure 3.8(a), for which kinematic interface compatibility conditions require that: (3.16) where (n = 1,2) represents a subset of the total DOFs that are in contact at a particular position. The compatibility condition that extracts the interface DOFs is described by a displacement operator of the form of a Boolean matrix represented as: 37 [ ]{ } (3.17) If however, a generalized bottom-up substructural synthesis approach is sought which may involve synthesis of substructures that are modeled separately with different mesh resolutions at interfaces, or by using dissimilar element types ? it would result in non-conforming mesh distributions at the contacting interfaces - as shown in Figure 3.8(b). The displacement operator will no longer have the form of a Boolean matrix in this case, and in order to enforce approximate geometric compatibility between substructures, Lagrange multipliers are generally introduced [46]. The Lagrange multipliers approximated by shape functions are not suitable for a position-dependent formulation, since, each time one substructure moves over another (say, by changing ? Figure 3.8), a different set of nodes come into contact at the interface, and employing the method in [46] would require a new set of Lagrange multipliers for each new position. To synthesize such substructures with incompatible meshes, an approximate model of surface interaction is obtained by defining a virtual condensation node placed at the center of each of the interface surfaces in contact, as shown in Figure 3.8(b), and subsequently enforcing displacement compatibility between these virtual nodes as: . (3.18) The displacement of the virtual node consisting of a set of translational and rotational DOFs ( is linked to the interface DOFs with a displacement operator such that: . (3.19) Coefficients of the displacement operator are extracted by linking the DOFs of the condensation node to the interface nodal DOFs it is meant to represent ( ) using a multipoint constraint (MPC) equation formulation. The MPCs represent the DOFs of the condensation node as a linear combination of all the DOFs of the nodes in contact at the interface [110]. 38 3.4.1. Multipoint Constraint Equation Formulation Flexible bodies can be loaded in a multitude of nodal DOFs simultaneously, making force interactions possible over the entire contacting interface. To avoid the issue of multiple solutions made possible by the condensation node representing only a subset of these simultaneously loaded DOFs, all interface nodal DOFs are represented by the condensation node DOFs using constraint equations without any interface reduction. The condensation node hence describes a net translational and rotational motion of the nodes of the interface surface it represents. All interface DOFs in contact are linked to the DOFs of the condensation node with the rigid MPC and/ or the interpolation MPC [111]. The formulation presented here, is for the generalized case of an individual substructural interface, hence the substructure superscript n is dropped. 3.4.1.1 Rigid Multipoint Constraint A rigid MPC formulation links the DOFs in contact such that part of the interface surface in contact and the condensation node move as a rigid system. For each interface node k being linked to the condensation node, the displacement is fully defined by the displacement and the orientation of the condensation node as [111]: (3.20) where is the vector from the condensation node to the node corresponding to the interface node k. The orientation of the node is represented by rotations around linear coordinate axes by assuming small angle approximations. The set of nodal DOFs for each node for volumetric elements in this work is represented as: . The coefficients of the displacement operator, , are extracted from the constraint equations setup using Eq. (3.20) and compatibility as in Eq. (3.18). The rigid MPC formulation introduces as many constraints as the interface DOFs it represents, less the DOFs of the condensation node. 3.4.1.2 Interpolation Multipoint Constraint Formulation The interpolation MPC formulation defines displacements and rotations of the condensation node as the weighted average of the motion of the interface nodes in contact. 39 The motion of the condensation node is fully described by the displacement of all interface nodes in contact for a total of p exterior DOFs as [111]: ? ? ? ? | | (3.21) where represents the weight factors for each DOF. To ensure that the condensation node represents the average motion of the contacting interface, the weight factors for each DOF are chosen proportional to the part of the interface surface its node represents; and are assigned as the coordinates of the nodes being coupled in this study. Since the modeling scheme is deliberately flexible to incorporate different substructures modeled using different mesh resolutions, the weights of all DOFs for an irregular mesh may not equal and may result in some DOFs being more heavily weighted than others. Moreover, since the moving body has its own frame of reference attached to it which is different than the global frame of reference, it may result in some DOFs having zero weights, which may result in numerical errors during condensation and synthesis. Hence, care should be exercised that each DOF being coupled has a real and non-zero weight associated with it. The coefficients of the displacement operator, , as in the case for the rigid constraint formulation are extracted from the constraint equations setup using Eq. (3.21) and compatibility as in Eq. (3.18). An interpolation MPC formulation introduces as many additional constraints as there are condensation DOFs. 3.5. Numerical Methods to Handle Constraint Equations On account of constraining the displacements at the interface, there is an increase in the potential energy of the system. In order to solve the constrained set of equations of motion, the potential function is modified by adding an extra ?energy? term and is solved by the Penalty method; and/ or the Lagrange multiplier method. 40 3.5.1. Penalty Method The variation of the modified potential functional in its standard form yields [106,112]: [ ] (3.22) where is the displacement operator whose coefficients are extracted from Eq. (3.20) or Eq. (3.21), and Eq. (3.18-3.19); is the penalty matrix, [ ] is a diagonal matrix of m penalty numbers corresponding to the m constraint equations, is the stiffness matrix, and is the external force vector. If in Eq. (3.22), the system of equations returns to the case of no constraints being imposed. As (i = 1,?,m) becomes very large, the penalty of violating constraints becomes large, so the constraints are closely satisfied. Numerical stability of this method is a function of the choice of , which are selected based on guidelines as described in [106,112]. If numerical stability is a concern over selection of proper penalty numbers, iterative and/ or optimal procedures for selection may be followed as suggested in [113,114]. 3.5.2. Lagrange Multiplier Method Using a discrete set of m Lagrange multipliers ( ) corresponding to the m constraint equations, and after taking the first variation of the modified potential functional we get [106,112]: [ ] { } { } (3.23) As opposed to the Penalty methods, the constraints are always satisfied in the Lagrange multiplier method. However, the size of the stiffness matrix increases by the factor of the number of Lagrange multipliers employed. 3.5.3. Substructural Synthesis The assembled undamped equation of motion ensuring compatibility at the interfaces for the two substructures in Figure 3.8 is represented as: [ ] { ? ? } [ ] { } { } (3.24) for the Penalty method, and 41 [ ]{ ? ? ? } [ ]{ } { } (3.25) for the Lagrange multiplier method; where { } represent the reduced substructural mass and stiffness matrices; and ( ) is the displacement operator coupling the substructures. Employing the floating frame of reference technique [111], these sets of equations, i.e. the eigenvalue problem form of Eq. (3.24) and (3.25), are solved for any desired position, by varying the position of the substructure 2, set by , as shown in Figure 3.9; thus obtaining the dynamic response at any position. Figure 3.9 Overall flow chart showing a sequence of operations to obtain substructurally synthesized position-dependent response 42 The proposed formulation presents a complete description of modeling the position-dependent dynamics based on defining constraint equations to synthesize reduced substructures. 3.6. Application ? Modeling Three Axis Milling Machine The position-dependency of the dynamic stiffness at the TCP of a three axis vertical milling machine is primarily due to the relative motion of the spindle-spindle housing moving over the vertical column. Hence, as a first step, only the spindle (including the tool, tool-holder), spindle housing, and column substructures ? shown in Figure 3.10 are modeled, reduced independently, and combined subsequently using constraint formulations. Figure 3.10 Substructural assembly of the spindle-spindle housing substructure with the column substructure through constraint formulations The reduced substructural model consist of interface DOFs complemented by a set of component generalized coordinates, as listed in Table 3.1. The interface DOFs for the 43 column substructure represent the guide-way top surfaces and surfaces in contact with the machine base; and, for the spindle-spindle housing substructure the interface DOFs represent the guide-block surfaces in contact with the column guide-ways and the DOFs corresponding to the spindle assembly including the tool-tool holder DOFs. The generalized coordinates correspond to the significant modes identified with the mode indicator function from the FRF simulated up to 10 kHz between the optimal driving and response locations. Uniform damping of the level of is assumed for all modes. In this way, 31 significant modes are identified for the spindle-spindle housing substructure, and 43 for the column substructure. Table 3.1 Division of DOFs for the substructural components Column Spindle-Spindle Housing Total Full Order Model 10908 15117 26025 Reduced Model Interface DOFs ( ) 762 1845 2681 Significant component modes ( ) 43 31 Table 3.2 compares the first 15 significant modes identified for the representative column substructure using the mode indicator function based sorting scheme with that of the standard CMS modal reduction scheme, i.e. in which the first few ( non rigid-body low frequency modes are retained. As is evident in representative results in Table 3.2, the mode sets are disjoint; and, that the mode indicator function based ranking scheme spans a wider frequency range than the standard modal reduction. For the standard modal reduction scheme to span the same frequency range as that of the mode indicator function based scheme, 310 mode sets would need to be retained for the column substructure, instead of the 43 presently identified as significant. In the case of the spindle-spindle housing substructure, 139 modes would need to be retained instead of the 31 identified as significant. The reduced set hence includes effects of higher order dynamics while keeping the reduced size to a minimum while automatically deciding the order of the reduced model. 44 Table 3.2 Top 15 mode subsets for standard modal reduction (SMR) based method and mode indicator function (MIF) based sorting for the column substructure, where is the natural frequency SMR Mode # 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 [kHz] 0.21 0.48 0.6 0.68 0.87 0.92 0.93 0.96 1.02 1.14 1.17 1.21 1.33 1.37 1.39 MIF Mode # 1 4 5 7 9 11 12 13 14 16 22 25 27 28 32 [kHz] 0.21 0.68 0.87 0.93 1.02 1.17 1.21 1.33 1.37 1.47 1.64 1.74 1.83 1.89 2.07 45 3.6.1. Verification of the New Improved Reduced Order Models To determine the accuracy of the improved reduced order model over the standard CMS scheme; and to verify the proposed improved variant; two criteria are employed: the normalized relative frequency difference; and, the modal assurance criterion. The normalized relative frequency difference (NRFD) is defined as: | | (3.26) where and are the eigenfrequencies of the full order model and the reduced order model respectively. The modal assurance criterion (MAC) in turn is defined as [108]: | | ( ) (3.27) where and are the mode shape vectors for the full and the reduced order model respectively. Generally, MAC values in excess of 0.9 should be attained for well-correlated modes and a value of less than 0.1 for uncorrelated modes. Both, the column and the spindle-spindle housing substructures were verified and the results are shown in Figures. 3.11-3.12. The figures compare the NRFD and MAC of the first 80 modes between the standard Craig-Bampton CMS method including standard modal reduction based sorting scheme, and that of the iterated improved component mode synthesis with mode indicator based sorting scheme. 46 Figure 3.11 Normalized relative frequency difference (NRFD) and modal assurance criterion (MAC) comparison for standard component mode synthesis scheme (left) and iterated improved component mode synthesis scheme (right) ? for the column substructure The reduced model response characteristics are more influenced by the choice and number of exterior/ interface DOFs than the geometry of the substructure. The more the number of interface DOFs, better approximated is the full model. For the example under consideration, the column substructure has less interface DOFs than the spindle-spindle housing substructure, hence making it more with diffficult to have approximated full order response. 47 Figure 3.12 Normalized relative frequency difference (NRFD) and modal assurance criterion (MAC) comparison for standard component mode synthesis scheme (left) and iterated improved component mode synthesis scheme (right) ? for the spindle-spindle housing substructure For the column substructure (Figure 3.11), the improved reduced order model has a frequency error of less than 1% for entire frequency spectrum considered, while the standard CMS method has errors as high as 60% for the high frequency modes. Similarly, the improved reduced order model is far better correlated to the full order model than the standard CMS scheme, having diagonal dominated near uniform unity MAC values; hence, making it the method employed for all subsequent analysis. The MAC for the improved reduced spindle-spindle housing substructure (Figure 3.12) has fewer off-diagonal terms as compared to the coulmn substructure (Figure 3.11), since its reduced set contains more exterior DOFs than the column substructure. 48 3.6.2. Structural Assembly for Different Positions of the Substructures Each of the four pairs of contacting interfaces between the four spindle-housing guide-block sets and the column guide-ways is represented by a pair of condensation nodes linked to the interface DOFs using multipoint constraint formulations of Eq. (3.20), and Eq. (3.21). The weight factors in Eq. (3.21) are assigned as the nodal coordinates. The modified (reduced) form of the synthesized constrained equation sets of Eq. (3.24) and Eq. (3.25) are solved using the Penalty method from Eq. (3.22) as well as the Lagrange Multiplier method from Eq. (3.23). To evaluate the reduced order substructurally synthesized position-dependent model, TCP FRFs simulated assuming uniform damping of the level of for all modes are compared with that of a full order model obtained from ANSYS? [104], wherein the substructures were ?glued? together. To evaluate the effect of the type of constraint formulation, results with different MPC formulations are compared in Figure 3.13. The rigid MPC formulation overestimates stiffness, which in turn results in higher eigenfrequency estimates [111]. The FRF for the rigid MPC seems to have shifted to right when compared with the full order model FRF, while the FRF with the interpolation MPC formulation reasonably approximates full order response, and hence is recommended as the method of choice to be employed for subsequent synthesis. Figure 3.13 Comparison of TCP FRFs for solution with different constraint formulations; assumed uniform damping of for all modes 49 Furthermore, results with different numerical solution schemes to the constrained problem, i.e. with the Lagrange multiplier method and the Penalty method are also compared in Figure 3.14. For penalty numbers selected as [ ] [112] both methods give the same response. Figure 3.14 Comparison of TCP FRFs for solution with different numerical methods; assumed uniform damping of for all modes As a second step, position-dependent response comparisons between the substructurally synthesized reduced order models and full order models are shown in Figure 3.15. The substructural synthesis being based on a floating frame of reference system, the moving frame attached to the TCP is adjusted relative to the reference frame by varying (in Figure 3.10). For each new position, a new set of constraint equations are established to describe the relations between the new nodes that have come into contact, by suitably updating the displacement operator, , within Eq. (3.24) or Eq. (3.25). Subsequently, the new reduced synthesized constrained equation sets of Eq. (3.24) or Eq. (3.25) are solved to obtain the position-dependent response. For the overall sequence of operations to obtain position-dependent response, see Figure 3.9 again. Three different positions are compared: the configuration when the spindle-spindle housing is at the top position, as shown in Figure 3.10; a mid position; and, a bottom position when the spindle-spindle housing combine has moved in the Z-direction by an amount of -0.2 m; and -0.4 m respectively. 50 Figure 3.15 Comparison of full order model and reduced order model TCP FRFs at three different positions: top position (top), mid position (middle), and bottom position (bottom); assumed uniform damping of for all modes As evident from comparisons in Figure 3.15, the substructurally synthesized reduced machine model is able to capture full order model results with little to no errors. Slight differences observed between the two models may be attributed more to the nature of the contacting surface interaction approximations using the constraint formulations than to the reduction process. The TCP FRFs are also observed to be strongly position dependent, with the lower frequency structural modes exhibiting stronger position-dependent behavior as compared to the higher frequency local tool-tool holder modes. 51 The position-varying dynamic behavior (natural frequencies and dynamic stiffness) is also tabulated in Table 3.3. As is evident, the global modes corresponding to the column (up to 150 Hz) and the spindle housing (up to 600 Hz) exhibit strong position-dependent behavior, varying by 6-10% in natural frequencies and by 50-75% in dynamic stiffness. The spindle-tool-tool-holder modes (600 ? 1200 Hz) on the other hand are more local in nature and do not exhibit strong position-dependency. The spindle mode around ~640 Hz varies by at most 4% in natural frequency and by ~15% in dynamic stiffness. The tool-tool-holder modes ~1000 - 1200 Hz have a negligible variation in natural frequency (< 1%) and at most a 12% change in dynamic stiffness. Table 3.3 Comparison of position-dependent behavior for the improved reduced model, where is the natural frequency, is the dynamic stiffness, Mode # Top Position Mid Position Bottom Position Associated mode shape [Hz] [N/?m] [Hz] [N/?m] [Hz] [N/?m] 1 73 2.8 77 3.7 81 4.8 Global column bending 2 124 4.4 125 4.1 128 3.8 Global spindle-housing + column bending 3 258 10.7 263 11.1 260 12.2 Global spindle-housing bending 4 310 11.6 318 12.6 325 12.5 Global spindle-housing bending + torsion 5 529 5.2 550 4.6 561 4.1 Global spindle-housing torsion 6 639 6.5 640 6.2 612 7.6 Spindle bending 7 853 12.1 - - 889 12.2 Spindle + tool-tool holder 52 bending 8 1066 4.1 - - 1083 3.6 Tool-tool holder bending 9 1149 2.1 1131 1.7 1148 2 Tool-tool holder bending The synthesized reduced model is able to reasonably approximate full model behavior and being ~1/10th the size of the full model (see Table 3.1) leads to considerable simulation time savings; taking ~10 seconds/ position as compared to ~6 hours/ position for the ?full? model (Intel? i3-380M processor with 4 GB RAM) thereby facilitating further position-dependent stability analysis. 3.7. Summary A systematic and computationally efficient procedure is proposed to model and evaluate the position-dependent dynamic behavior of a three axis milling machine tool based on substructural synthesis of improved reduced order models as an alternative to presently used time consuming full FE models by the designers. Machine substructural components were reduced and subsequently synthesized using adaptations of constraint equations which tolerate mesh incompatibility at their contacting interfaces; allowing for modular design. Substructural components were reduced using an improved variant of the component mode synthesis method to retain interface DOFs complemented by a set of judiciously selected component modes using a novel mode selection criterion. The mode sets thus identified are able to represent higher order dynamics of the substructure while keeping the order of the reduced model to a minimum by spanning a much wider frequency range with fewer modes than would be required with standard CMS methods. Substructural syntheses with the interpolation MPC formulation were found to be better than the rigid MPC formulation. For a correct choice of the penalty number, both numerical solution techniques to the constrained problem provided the same results. Dynamic response characteristics of the substructurally synthesized reduced order FE models were verified against corresponding full order models. 53 The substructures were rigidly coupled in the preliminary investigation, without having modeled the stiffness and damping at the joints. Joint characteristics severely affect tool point response ? and are treated separately in Chapter 4. Joints once modeled, the reduced model substructural synthesis formulation developed here is extended to model a complete machine tool in Chapter 5. 54 Chapter 4 Modeling and Identification of Machine Tool Joints1 4.1. Overview The multitudes of joints in a typical machine tool have varying degrees of influence on the response of the machine at the point of interest, i.e. at the TCP. The low frequency TCP behavior is severely affected by the joints between the various structural substructures such as the joints between the column and the base, and the base and ground; whereas the higher frequency TCP behavior is more dependent on the joints between the tool and tool holder as well as between the tool holder and the spindle. To obtain the full frequency spectrum response a two pronged strategy is proposed in this chapter ? as shown in Figure 4.1. At first, each of the major machine tool components for the machine under consideration, i.e. a three axis vertical machining center ? FADAL 2216 are all modeled independently as discussed in Chapter 3. These major substructures (base, column, table, and spindle-spindle-housing) are synthesized using a simplified approach of idealizing the contacting interfaces between these substructures as being connected by linear spring elements ? as discussed in Section 4.2. The response at the spindle nose obtained from the synthesized machine model (Substructure III ? Figure 4.1) is further combined with the tool-tool-holder combine (Substructures I + II, Figure 4.1). Since the response at the spindle nose is available from the virtual machine tool model, measurement of rotational and translational FRFs at the spindle nose, which can be challenging is no longer necessary. Moreover, the virtual machine tool model enables us to include the effects of machine tool structural dynamics on the identified joint dynamics. As a second step, the tool-tool-holder assembly response is obtained by synthesizing the response of the tool (Substructure I) with the tool-holder response (Substructure II) using the IRC approach ? as discussed in Section 4.5. The IRC method can be used to identify interface characteristics between the tool and the tool-holder and between the tool-holder and the spindle, since these substructures are easy enough to disassemble and measure in their 1 Parts of this Chapter have been published as Mehrpouya, M., Law, M., Park, C., Park, S. and Altintas, Y., 2013, Joint Dynamics Identification on a Vertical CNC Machine, 2nd International Conference on Virtual Machining Process Technology, Hamilton, Canada, [2] 55 free-free configuration. This two-stage substructural synthesis approach allows different substructural combinations to be modeled independently and coupled subsequently. Figure 4.1 Two-stage substructural synthesis of the machine tool. Tool (Substructure I) is synthesized with the tool-holder response (Substructure II); and, the synthesized tool-tool-holder combine is subsequently coupled to the response of the synthesized substructural machine model (Substructure III). 4.2. Modeling Joints between Structural Substructures Joint characteristics for bolted connections between structural members; connections between the guide-block and rail and between ball-screw and nut; and, bearings; depend on 56 parameters like preloads, contact surface conditions, bearing types, friction, and damping. The information about joint stiffness and damping is seldom available at the design stage; making joint modeling non-trivial in the FE environment. Treatment, modeling and identification of joint characteristics are independent of the machine being represented by its full order FE model or its reduced order FE model. Hence, as a first step, joints are introduced between the structural components for the full order FE model in this chapter; and, are extended to the complete substructurally synthesized reduced machine model in Chapter 5. Joint flexibilities are idealized by spring elements. This approach, also referred to as the ?whole-joint? approximation [14] involves imposing simplified kinematics across the joint interface, relating the kinematics of all the nodes on each side of the interface to the kinematics of a representative node (either real or virtual), which is subsequently coupled to a representative node on the other contacting interface with an idealized spring. 4.2.1. Idealizing Bolted Connections Though there exist several fasteners between various substructures of the machine, only those connections that may contribute towards the overall tool point compliance are modeled in this study as connected by linear springs. As such, only the major interfaces between the base and column, between the spindle-housing and the column, and between the base and the ground are modeled; while other connections are assumed to be in rigid contact. Bolted connections are modeled as linear spring (bar) elements with stiffness as a function of material and geometric properties of the fasteners. The spring stiffness for bolted connections is expressed as: [ ] (4.1) where is the nominal cross-sectional area of the bolt; is the modulus of elasticity (steel assumed); and ? is the length of the bolt. 4.2.2. Idealizing Bearings as Springs Complex analytical models for bearings have been developed by several past researchers which depend on defining bearing stiffness as a function of applied force, resulting deformations, rotational speed, and even including centrifugal and gyroscopic effects 57 [115,116]. However, for the sake of simplicity, bearings when preloaded may be represented by a stiffness value that is usually available for various bearing types and arrangements from the respective manufacturers? catalogues. For each of the three ball-screw feed drive models in the machine being model, the ball-screws are supported with the help of angular-contact ball bearings with a medium-high preload setting of ~700 N. For these preload levels, the axial stiffness of the bearings is taken as 135 N/?m, and the radial stiffness is taken as 95 N/?m. For the case of the spindle bearings, for a medium-high preload setting of ~1000 N on the angular-contact ball bearings, the radial stiffness is taken as 212 N/?m, and the axial stiffness is taken as 97 N/?m. Additional details on modelling the support bearings are described in [27,40]. 4.2.3. Idealizing Connections between Guide-ways and Guide-blocks, and between Ball-screw and Nut The interfaces between the guide-ways and guide-blocks and between the ball-screw and nut are idealized as being connected by linear spring elements using the ?whole-joint? approximation with the equivalent contact stiffness values obtained from manufacturers? catalogues. Joints at these contacting interfaces are idealized as two translational springs perpendicular to the direction of motion with no resistance, i.e. no spring in the direction of motion. For the machine under consideration, each axis has two guide-ways and four guide-blocks; and one ball-screw nut interface respectively. Equivalent contact stiffness for each of the three axes for the guide-block and guide-rail interface is assigned as 187 N/?m (THK SVR series); and as 280 N/?m for the ball-screw-nut interface (THK SBN series) [117]. 4.3. Validating Model Predicted Response at Spindle Nose with Measured Response At first, model response at the spindle nose obtained by modeling the joints between the structural substructures as described above is compared with response obtained with the case of all the structural components assumed to be in rigid contact by assuming uniform damping for all modes. Damping is subsequently measured from experiments and is used to update 58 model predicted response, following which updated model predicted response is compared with measured behavior. 4.3.1. Comparison of Model with Rigid and Spring Connections Spindle nose FRFs in machine principal directions, i.e. in the X and Y directions are compared in Figure 4.2, for the case of model with rigid and spring connections. Response for the full order models is obtained by carrying out modal analysis for the full FE model within the ANSYS? environment. Comparisons are limited to the low-frequency regime of up to 400 Hz, i.e. to the frequency range influenced by the machine design dependent structural components but not as application dependent tools and tool holders. Figure 4.2 Response comparisons at spindle nose for machine with rigid connections and spring connections in X (top) and Y (bottom) directions; assumed uniform damping of for all modes 59 As evident from comparisons in Figure 4.2, the response at the spindle nose for the case of rigid connections results in higher frequency estimates than the response with spring connections modeled as described in Section 4.2. Also, the dynamic stiffness for the case of connections with springs is lower than the case of rigid connections. The main influence of joints is to bring about a reduction in the stiffness of assembled structure(s), which in turn has the effect of lowering the natural frequencies and the dynamics stiffness, as is observed in Figure 4.2. A uniform damping ratio of ? = 0.05 has been assumed in simulating the above FRFs. These damping levels are updated by correlating measured response with model predicted response - as presented in the following section. 4.3.2. Experimental Modal Analysis The test setup to identify modal damping on the machine being modeled is shown in Figure 4.3. An instrumented impact hammer is used to excite the machine at the spindle nose and a laser vibrometer is used to measure the displacement response. The measured FRFs are curve fit within CUTPRO? [118] to identify modal damping levels. The first few dominant modes of the machine that are measured are listed in Table 4.1, along with the identified modal parameters. The identified modal parameters (modal damping levels) are subsequently used to update the model predicted response. Figure 4.3 Experimental test setup for modal analysis on FADAL 2216 ? three axis vertical machining centre 60 Table 4.1 Modal Parameters obtained from measurements on FADAL 2216, where is the natural frequency, is the dynamic stiffness, and ? is the damping ratio X direction Mode # [Hz] [N/?m] ? [%] 1 36 23.1 6.4 2 97 12.5 6 3 130 16.9 4.5 Y direction 1 26 7.5 6 4.3.3. Comparison of Measured and Updated Model Response Measured spindle nose FRFs computed with the modal parameters listed in Table 4.1 are compared with updated model predicted response in Figure 4.4, for both the machine X and Y directions respectively. Measured and model predicted behavior of the dominant modes is also compared in Table 4.2, which also lists the errors in prediction of the dynamic stiffness?s and natural frequencies rounded to the nearest integer. Figure 4.4 Measured response comparisons at spindle nose with model predicted response for rigid connections and spring connections in X (top) and Y (bottom) directions. 61 Table 4.2 Comparison of dominant low frequency modes for model predicted response (with/ without joint effects) with measured response at the spindle nose, where is the natural frequency, is the dynamic stiffness, is the error in natural frequency prediction and is the error in dynamic stiffness prediction Measured Model with joints Rigid model X Mode # [Hz] [N/?m] [Hz] [%] [N/?m] [%] [Hz] [%] [N/?m] [%] 1 36 23.1 32 -12 8.7 -165 42 14 15 -55 2 97 12.5 68 -42 4 -210 81 -20 5.7 -120 Y 1 26 7.5 31 16 10 25 56 55 29 72 As is evident from Figure 4.4 and Table 4.2, even though model predicted response with/without joints does not have the same number of modes as the measured model, model predicted response is able to capture the general trend of the measured response. Though model with joints is able to quite accurately capture the Y directional behavior, both in frequency matching and dynamic stiffness approximations as compared to the model with rigid connections; there exist considerably errors in approximating the X directional response; particularly in the second dominant mode in the X direction. In the respect of frequency matching alone, model with joints is better able to approximate measured modes with errors in predicting frequencies ranging from 12% to 16% at most, with the exception of the second dominant mode in the X direction, for which the error in approximation is as much as ~42%. Model with rigid connections have a higher error in predicting frequencies ranging from 14% to 55%. These errors in the low-frequency modes, as well as the order mismatch between the model predicted behavior and measured behavior may partially be attributed to the modeling simplifications in representing some minor machine components like the automatic tool changer and cabinets by lumped mass elements. Moreover, since the true nature of the contacts between these substructural interfaces is more complex than can be captured with an idealized spring connection, it may further explain the high level of errors in approximating the dynamic stiffness, ranging from 50% to 120% for the model with rigid connections, and slightly higher error for the model with 62 joints, ranging from 25% to 210%. These higher error levels in dynamic stiffness estimation for the model with joints are carried over from the differences observed in Figure 4.2, i.e. due to the inclusion of joints in the model that result in lowering the overall system stiffness. These errors may also be attributed to the difficulty in approximation of the base mounting pads by linear springs whose contact stiffness if underestimated, lead to errors in accurate prediction of the dynamic stiffness of the low-frequency modes. These errors may be minimized by modeling the machine and its joints with higher degree of fidelity. Overall, even though the model predicted response still has considerable errors, the model with joints better approximates the measured behavior as compared to the rigid model. Since these models are meant to aid in engineering decision making, as a first level of approximation for subsequent analyses and investigations, the error levels though high, are deemed acceptable. 4.4. Finite Element Models of the Tool and Tool-holders The tool and tool-holder response is to be synthesized with the response at the spindle nose obtained above. Since the spindle assembly was modeled with Timoshenko beam models as described in Chapter 3, to ensure element type compatibility, the tools, and tool-holders are also modeled with Timoshenko beam models. Each of the FE models of the tools and tool-holders are all checked for convergence and are subsequently validated with their physical counterparts by comparing response characteristics (FRFs) in their unsupported (free-free) configurations. Material properties assigned for the tools that are made of carbide are: modulus of Elasticity of 550 GPa; and, density of 15630 kg/m3. For the tool-holder made of steel, material properties assigned are: modulus of Elasticity of 210 GPa; and, density of 7800 kg/m3. Since the joint identification procedure discussed below also involves use of cylinders (representing blank tools) to identify and validate the joint characteristics, the cylinders that are made of steel are modeled with Timoshenko beam elements. The response at the spindle nose obtained from the validated virtual machine tool model is combined with the tool-tool-holder response with the inverse receptance coupling (IRC) approach to obtain the TCP FRF ? as discussed in the next Section. To avoid measuring rotational receptances which was very challenging and difficult, only translational FRFs of 63 the assembled structures are used in the identification. For the substructures, both rotational and translational FRFs are extracted from the FE models and used in the identification. 4.5. Identifying Joints between Substructures First, the finite element (FE) model of the tool and the tool-holder along with the measured frequency response functions (FRFs) of the free-free tool-tool-holder assembly are used to find the joint?s dynamic properties between the tool and the tool-holder. Secondly, the validated FE model of the machine tool is employed in the IRC method to obtain the joint?s dynamic properties between the holder and the spindle using the measured receptances on the actual physical machine tool structure. The identification of joint dynamics between the tool-holder and the spindle is performed by using the measured receptances along the tool while the tool-tool-holder assembly is inserted inside the spindle. 4.5.1. Inverse Receptance Coupling Let the substructures A and B, shown in Figure 4.5, be connected through a joint which is comprised of rotational and translational elements. These two substructures may represent any of the three substructures in Figure 4.1. Points nA and nB represent the internal locations on Substructures A and B, and, points cA and cB represent the DOFs connected through the joint section. Figure 4.5 Substructures coupled with a joint The relation between the displacements and the forces in each substructure can be defined as: 64 S SS SS S SJS S SJS S S Snn nn nc ncn nS S S Sn nn nn nc nc nS S S Sc ccn cn cc ccS S S Sc ccn cn cc cch l h lx fn o n o Mx f fh l h lM Mn o n o??? ?? ? ? ?? ?? ? ? ?? ?? ? ? ??? ? ? ?? ??? ? ? ?? ?? ? ? ?? ??? ? ? ?? ? (4.2) where xiS and ?iS, (S = A, B; i = n, c), represent the translational and rotational displacement vectors at location i on substructures A and B; and, FiS = {fiS , MiS} represents the vector of force and moment at location i on each substructure. The receptance components are defined as hij= xi/fj, lij= xi/Mj, nij= ?i/fj and oij= ?i/Mj; where represents the displacement-to-force receptance; represents the displacement-to-couple receptance; represents the rotations-to-force receptance; and represents the rotation-to-couple receptance. The equilibrium condition at the joint part is: 0??????????????JBJBJAJAMfMf (4.3) where FSJ = {fSJ , MSJ} is the vector of force and moment in the joint section at the connecting locations to substructure S. Using the equilibrium conditions, the equation of motion at the joint part can be written as: ? ???????????????????????????????????????????????? ?JAJAJJAJAxxJAJAJrrJttBcABcAMfHMfickickMfhhxxCc10000???? ?? (4.4) where HJ denotes the receptance matrix of the joint and subscripts t and r represent the translational and rotational directions respectively. Considering that forces and displacements of the internal coordinates do not change before and after coupling, i.e. xnA=xn; ?nA=?n; xnB=xn; ?nB=?n; FnA=Fn; MnA=Mn; FnB=Fn; and, MnB=Mn; the assembled structure?s FRFs can be obtained by inserting Eq. (4.4) into Eq. (4.2) as [2,93]: ?????????????????????????????????????????nBcBcAnAnBnBnBcBnBcAnBnAcBnBcBnBcBcAcBnAcAnBcAcBcAcAcAnAnAnBnAcBnAcAnAnAnBcBcAnAFFFFGGGGGGGGGGGGGGGGXXXX (4.5) 65 where [ ] represents the assembled structure?s FRFs, and Xi = {xi ?i}T is the displacement vector of the assembled structure. Considering that two internal locations are available for the measurements on a structure, nA = {1,2}, and naming the connecting points on substructures A and B as 3 and 4, two assembled structure?s FRFs are expanded as: 111, 11, 13 13 33 33 31 3111 11 44 4411, 11, 11 11 13 13 33 33 44 44 31 3112, 12, 12 1212, 12, 12 1200Jtt tr ttJrt rr rrtt trrt rrG G h l h l h lhh l h lG G n o n o n o n o n ohG G h lG G n o?? ?? ?? ? ? ? ? ? ? ?? ? ? ?? ? ? ?? ?? ?? ? ? ? ? ? ? ?? ? ? ?? ?? ?? ? ? ?? ? ? ? ? ?? ? ? ?? ?? ? ? ??? ? ?? ?? ?113 13 33 33 32 3244 4413 13 33 33 44 44 32 3200JttJrrh l h l h lhh ln o n o n o n oh?? ?? ?? ? ? ? ? ?? ?? ? ?? ?? ?? ? ? ? ? ?? ? ?? ?? ?? ?? ? ? ? ? ?? ?? ? (4.6) Since the only assembled structure?s FRFs that can be measured easily are translational FRFs, G11,tt and G12,tt are expanded as: 111 12 12 21 12 12 211112 12 12 22 12 12 2221[( ) ( ) ]11,1[( ) ( ) ]12,xG h h b l b h h b l b nrr rt tr tttt f b b b brrtt rt trxG h h b l b h h b l b nrr rt tr tttt f b b b brrtt rt tr? ? ? ? ? ? ??? ? ? ? ? ? ?? (4.7) where 33 33 44 4433 33 44 4400Jtt tr ttJrt rr rrb b h l hh lb b n o n o h? ?? ? ? ? ? ?? ? ? ? ?? ? ? ? ? ? ? ?? ?? ? ? ? ? ? (4.8) The two explicit joint parameters of httJ and hrrJ are obtained symbolically by simultaneously solving these two equations in Eq. (4.7) using MATLAB? symbolic toolbox, to yield: ? ?4433231123,11131213,1232213313213321132,11311231,12 hhlhlGlhlGnhnhhhhbhGbhhbhGbhttttrtttrtrtttrtJtt ??????????? ? ?12, 31 12 31 32 13 31 11, 32 11 32 31 13 3233 4412, 31 12 31 11, 32 11 32tr tt tr tr tt trJrrtt ttb G n b h n h l n b G n b h n h l nh o oG h h h G h h h? ? ? ? ?? ? ?? ? (4.9) Based on Eq. (4.9), the joint FRFs are obtained by two assembled structure?s receptances and the substructures? FRFs which are obtained through the FE models. This eliminates the 66 dependence of identification on the rotational FRF measurements and gives an explicit solution for the joint FRFs. Equation (4.9) will be repeatedly used in the identification of joint dynamics in the following sections. At each identification stage, two measured receptances along with the validated FE models of substructures are used to find the joint?s FRFs. In order to assess the accuracy of identified joint?s FRFs, different tools are inserted in the tool-holder and the identified joint?s FRFs are used in Eq. (4.7) to build the new structure?s FRFs. The reconstructed FRFs are then compared with the directly measured FRFs on the new structure. 4.5.2. Joint Identification between the Tool and Tool-holder The IRC procedure necessarily requires measurements on the substructures in their assembled configuration. The experimental setup to obtain these measurements is shown in Figure 4.6, and includes a CAT40 tool-holder, two cylinders, 50 mm and 70 mm long, and one actual end mill, 90 mm long, which were each inserted 30 mm inside the tool-holder, as also schematically represented in Figure 4.7. The joints are identified for a given tool and tool-holder combination, and subsequently validated with different tools. The 70 mm long cylinder is used in the identification step, and the 50 mm long cylinder along with the 90 mm long tool is used in the validation step. These components, i.e. the tool-holder and each of the cylinders and the actual tool were all modeled with Timoshenko beam elements in the FE environment. The translational as well as rotational FRFs at locations H44 and H45 on the tool-holder (see Figures 4.6-4.7) were obtained from its FE model. The cylinder/ tool FRFs, H11, H12, H22 and H13 were also obtained from their corresponding FE models. Substructural FE models for each of the tool(s) and the tool-holder were all validated against measurements in free-free configurations, i.e. unsupported conditions. Damping levels for the individual free-free FRFs for the FE models are obtained from measured modal damping, by assuming a Rayleigh definition of proportional damping [108]. 67 Figure 4.6 Free-free test setup for the tool-tool-holder combination Figure 4.7 Schematic of different tool and tool-holder assemblies The procedure to identify joint dynamics between the tool and the tool-holder is as shown in the flow chart in Figure 4.8. The flow chart depicts the identification procedure which was done on the 70 mm cylinder and the validation step which was done on the 90 mm tool and the 50 mm cylinder. For identification of the joint between the shank of 70 mm cylinder and the tool-holder, two measurements were performed at locations 1 and 2 (see Figures 4.6-4.7) with the tool-tool-holder combine in free-free conditions. This information, along with the 68 FRFs of the tool-holder and the cylinder, i.e., H44, H11, H12, H13, H22, H23 and H33 all obtained from the respective FE models are inserted into Eq. (4.9) to find the translational and rotational FRFs of the joint. The identified joint FRFs are shown in Figure 4.9 ? for the translational FRF, httJ, and in Figure 4.10 ? for the rotational FRF, hrrJ,. Figure 4.8 Flow chart showing procedure for joint identification, and validation between the tool and the holder 69 Figure 4.9 Identified translational joint FRF, httJ, between the tool and the tool-holder Figure 4.10 Identified rotational joint FRF, hrrJ, between the tool and the tool-holder From these two figures, i.e. Figures 4.9-4.10, the structural modes of the joint are observed to be between 5 kHz to 6 kHz. If the dynamics of a structure which uses this tool-tool-holder setup is sought around these frequencies, i.e. between 1 - 10 kHz, the effects of the joint between the tool and the tool-holder cannot be ignored. To investigate the accuracy of the identified joint FRFs, and, the potential improvements in the assembled structure FRFs resulting from the joint FRFs, the FRFs for the 50 mm cylinder and 90 mm tool were reconstructed by inserting the identified joint?s FRFs and substructures? FRFs in Eq. (4.7). The reconstructed free-free assembly FRFs are compared with the measured FRFs along with the case for the substructures being rigidly coupled in Figures 4.11 and 4.12. 70 Figure 4.11 Direct FRFs for the tool-tool-holder with 50 mm cylinder (tool), G11_50mm. Figure 4.12 Direct FRFs with tool-tool-holder with 90 mm tool, G11_90mm As evident in Figures 4.11 and 4.12, a considerable improvement in the prediction of the tool-tool-holder assembly was obtained with identified joint dynamics compared to the rigid joint connection. The error in prediction for first two natural frequencies has improved from 24.0% and 9.7% in the rigid joint assumption to 5.0% and 2.1% in the reconstructed FRF. The high level of noise observed in the reconstructed FRFs is due to the use of the measured response during the identification step. Though a Savitzky?Golay filter [119] was applied to the recorded signals before analysis to filter out some of the measurement noise, any residuals leftover may magnify the error while dealing with the matrix inversion in the 71 identification step, in Eqs. (4.7) and (4.9). This may also explain the additional modes observed in the reconstructed FRFs in Figures 4.11 - 4.12. At this tool-tool-holder assembly level, the substructures are symmetric, hence the X and Y directional response is the same. When this holder-tool assembly is inserted inside the spindle, the interface between these two bodies influences the dynamics at the TCP. The joint at this interface is identified in the next section. 4.5.3. Joint Identification between Spindle and the Tool-holder Two substructures are considered in this section: a FADAL 2216 machine tool whose response at the spindle nose was obtained from the validated virtual machine model in Section 4.3; and, the tool-tool-holder combine whose mathematical model was built by considering the joint effects between tool and the tool holder ? as discussed in Section 4.5.2. Having a validated virtual model of the machine center reduces the need to perform several measurements at the spindle nose that was required in earlier studies to obtain translational and rotational FRFs [93,96,120]. The joint identification procedure was done for the case of the 70 mm cylinder and tool-holder assembly clamped at the spindle, as shown schematically in Figure 4.13. Joint dynamics once identified, TCP FRFs for a tool with 90 mm overhang length from the tool holder collet face was constructed and validated against measurements ? using the procedure outlined in the flow chart in Figure 4.14. Figure 4.13 Schematic of spindle, holder, tool assemblies 72 Figure 4.14 Flow chart showing procedure for joint identification, and validation between the spindle nose and the tool-tool-holder combine Two measurements were made along the cylinder ? shown schematically in Figure 4.13; one at the cylinder tip and another 20 mm away from the tip, to obtain the G11,tt and G12,tt FRFs. These FRFs along with the FRFs at the spindle nose, H44, (see Figure 4.13), and the FRFs for the holder-cylinder assembly, H11, H12, H13, H23 and H33, were inserted into Eq. (4.9), shown in the flow chart in Figure 4.14 to obtain the joint?s translational and rotational FRFs, httJ and hrrJ. Figure 4.15 shows the translational joint?s identified FRF, and, as evident several structural modes exist in the joint FRF. This shows that the joint between the spindle and the tool-tool-holder combine has more significant effects on the dynamics of the assembled structure than the joint between the tool and the tool-holder, which showed only one structural mode (Figures 4.9-4.10). 73 Figure 4.15 Identified translational Joint FRF between Spindle and tool-holder, httJ. To further validate the accuracy of the identified joint properties, the joint?s FRFs were used to construct TCP FRFs for the case of the tool-tool-holder combine with a tool of 90 mm length. To do this, the 70 mm cylinder-tool-holder assembly was replaced with the 90 mm tool-tool-holder assembly. The spindle nose FRFs, H44, the tool-holder FRFs, H11, H13, in Figure 4.13, and the joint?s FRFs, hJtt and hJrr, were inserted into Eq. (4.7) to find TCP FRFs. Figure 4.16 shows the comparison between the reconstructed FRF, the measured FRF at the tip of 90 mm tool and the assembled FRF obtained by considering a rigid joint between the tool-holder and spindle, hJtt = hJrr = 0 in Eq. (4.7). Figure 4.16 Direct FRF at TCP with 90 mm tool, G11_90mm. 74 The improvements that were obtained in predicting TCP FRFs compared to the rigid joint approximation showed the importance of the joint dynamics properties both between the tool-holder and the tool and between the tool-holder and the spindle. If an accurate prediction at the TCP is sought, the joint dynamics effects should be taken into account at both places. The proposed method in this study had several assumptions and limitations in the modeling, experiments and simulations. First, the off-diagonal terms in the FRF matrix are assumed to have negligible effects on the assembled structure FRFs. This assumption is true if, as in this study, the joint section mainly acts as a connecting element and imposes stiffness and damping to the structure. The behaviour of the joint was also considered to be linear in the studied frequency range, so the friction was modeled with a linear viscous damping element and nonlinear effects such as slipping were ignored [121]. The behaviour of the joint was also considered to be time-invariant and stable. From the modeling point of view, differences between the measured and simulated response at the spindle nose will compound inaccuracies in the identification. Idealization of cylinders and tools with Timoshenko beam models may have also caused some deviations in the identified joint properties. Lastly, the accuracy of the identified joint dynamics and their validity depend on similarity of joint conditions in the identification and validation structures. Many conditions such as contact area, pre-stress and manufacturing tolerances can affect the joint behaviour. The proposed methods are applicable when the influential conditions on the joint?s dynamic behaviour remain constant after replacing different substructures. 4.6. Summary A systematic procedure to model and identify joint characteristics was presented in this Chapter. A two-stage substructural synthesis approach was proposed to facilitate response analyses over the full frequency range of interest. The low frequency behavior was approximated by idealizing joints between the various structural components of the machine as linear springs, for which damping was incorporated by correlating model predicted response with measured response. Having a validated virtual model of the machine center reduces the need to perform several measurements at the spindle nose that was required in earlier studies to obtain translational and rotational FRFs. The mid-to-high frequency 75 behavior was approximated by applying an inverse receptance coupling method which identified joint dynamics between the tool and the tool holder and also between the tool-holder and the spindle. Though there still exist some errors between model predicted response and measured behavior, a considerable improvement in response prediction for the assembled structure was observed by modeling and identifying the joint dynamics as compared to treating the joint as rigid. These techniques of modeling and identification of joint characteristics are extended in the following Chapter for the case of modeling the position-dependent dynamics of a complete three axis vertical milling machine using the substructurally synthesized reduced order modeling scheme proposed in Chapter 3. 76 Chapter 5 Rapid Evaluation and Optimization of Serial Machine tools with Position-dependent Dynamics and Stability1 5.1. Overview To facilitate rapid evaluation and optimization of the dynamic behavior of machine tools, the reduced model substructural synthesis approach proposed in Chapter 3 is extended in this Chapter to model and assess the position-dependent dynamic behavior of a complete representative three-axis vertical machining center ? FADAL 2216. Effects of joints are also included by incorporating joint modeling and identification techniques presented in Chapter 4. The complete iterative virtual integrated position-dependent process-machine interaction approach to assess machine performance is shown in Figure 5.1. Figure 5.1 Flow chart of a virtual integrated position-dependent process-machine interaction approach for designing milling machines ensuring targeted productivity. 1 A version of this Chapter has been published as Law, M., Altintas, Y. and Phani, A. S., 2013, Rapid evaluation and optimization of machine tools with position-dependent stability, International Journal of Machine Tools and Manufacture, Vol. 68 [3]. 77 At first, the machine CAD model is partitioned into its major substructures. These models are modeled with FE and reduced with the improved reduction scheme proposed in Chapter 3. The reduced models are synthesized at the desired position to obtain position-dependent dynamic response. The generalized reduced model substructural synthesis formulation is discussed in Section 5.2. This model is verified by comparing with corresponding full order model results and validated against measurements; both of which are discussed in Sections 5.3 and 5.4. The machine under consideration is envisaged as a machine capable of high engagement heavy-duty steel machining (milling); hence, a performance requirement of a minimum of 3 mm chatter-free stable depth of cuts within the entire work volume is defined as the criterion against which to evaluate the machine tool. This criterion also identifies parameters limiting the target productivity levels. For a defined set of cutting conditions and productivity levels, the position and feed-direction-dependent stability is evaluated in Section 5.5. Mechanical parameters limiting productivity are identified and modified to meet target productivity, which forms the discussions in Section 5.6. 5.2. Generalized Substructural Formulation for Position-dependency Position-dependency at the TCP is modeled based on a two stage substructural assembly approach. At first, each of the major substructures of the machine under consideration, namely: spindle-spindle-housing, column, base, cross-slide, and table are modeled and reduced independently. These reduced substructural models are subsequently synthesized together with the three individual feed drive models as shown in Figure 5.2, allowing for efficient prediction of position-dependent response at the spindle nose. This is followed by a second stage substructural assembly procedure which involves coupling the tool-tool-holder response (Substructure I) to the response obtained at the spindle nose from the first stage (Substructure II) using a receptance coupling approach while including the joint dynamics as discussed in Chapter 4, Section 4.5. Joint characteristics between the tool and the tool-holder are identified as previously discussed in Chapter 4, Section 4.5.2; and hence are not discussed further here. This two-stage substructural synthesis facilitates modularity in the design and evaluation process by allowing different design variants of the structural components to be combined with different tool-tool-holder sub-assemblies. 78 Figure 5.2 Two stage substructural synthesis of the machine tool. Tool-tool-holder response (Substructure I - right) is coupled to the position-dependent response of the synthesized substructural machine model (Substructure II - left). 5.2.1. Substructure Model Reduction FE models of the structural substructures, each of the three independent feed-drive models, and the spindle assembly are all modeled as discussed in Chapter 3, Section 3.2. Each of these models is exported, after convergence checks, to the MATLAB environment for reduction and synthesis as discussed in Chapter 3, Sections 3.3-3.4. The size of each of the reduced models which consists of the exterior/ interface DOFs, , complemented by a set of P significant modes is listed in Table 5.1. 79 Table 5.1 Division of DOFs for individual substructures Spindle-Spindle Housing Column Base Cross-slide Table Total Full order model 21992 45983 97214 50201 15623 231013 Reduced Model Interface DOFs ( ) 1037 1533 1945 2598 1908 9173 Significant component modes (P) 27 28 46 27 24 5.2.2. Substructural Synthesis Synthesizing each substructural interface with a set of constraint equations, the undamped equation of motion for the synthesized reduced model at a particular position is: [ ] { ? ? ? ? ? ? } [ ] { } = { } (5.1) where { } are the reduced substructural mass and stiffness matrices and the subscripts and correspond to the spindle-spindle housing combine, column, base, cross-slide, and table respectively. ( ) is a displacement operator whose coefficients are obtained by ensuring displacement compatibility at the contacting interfaces with an interpolation constraint formulation as discussed in Chapter 3, Section 3.4.5. in Eq. (5.1) represents a column vector containing a discrete set of Lagrange multipliers corresponding to the number of constraint equations [106]. 80 The synthesized model enables prediction of dynamic response at the spindle nose as one component changes its position relative to another by solving the eigenvalue problem in Eq. (5.1), by varying the tool position in the work volume by adjusting or (see Figure 5.2). For each new position, since a new set of nodes come into contact, while others fall out of contact, the displacement operator in Eq. (5.1) is updated by instantaneously coupling/ de-coupling the corresponding nodes on the interfaces. Since the substructurally synthesized reduced machine model size is ~1/25th the size of the full model (Table 5.1), it allows for very efficient position-dependent dynamic modeling of the machine. Further modularity in the design process is facilitated by synthesizing the frequency response function (FRF) at the spindle nose with that of separately modeled tool-tool-holder model response as discussed in the next section. 5.2.3. Tool Point Response Predictions with Receptance Coupling The component receptances for the tool and the tool-holder that are modeled with Timoshenko beam elements can be expressed in a compact matrix generalized form as: (5.2) where is the generalized receptance matrix that describes both translational and rotational component behavior, and and are the respective measurement and excitation locations and, and are the corresponding generalized displacement/ rotation and force/ couple vectors. is made up of the component receptances as in Eq. (4.2) in Section 4.5.1. The direct receptances at, and cross receptances between the free-end (i.e. location in Figure 5.2) and the coupling end (location in Figure 5.2) of the tool-tool-holder model are: (5.3) (5.4) . (5.5) Due to symmetry and Maxwell?s reciprocity, . The direct receptances at the spindle nose (location 2 in Figure 5.2) are similarly represented as: (5.6) Ensuring equilibrium conditions at the joint part, i.e. at the connection between the tool-tool-holder and the spindle nose, the assembled receptances, , at the TCP in the generalized form are given as: 81 ( ) (5.7) where has the same structure as , and its constituent receptances are obtained from solutions in Eq. (5.3) ? (5.6). within Eq. (5.7) has the structure as in Eq. (4.9), as previously discussed in Section 4.5.1; and its joint characteristics are identified from measurements on the physical test machine in assembled configuration, i.e. with the tool-tool-holder inserted inside the spindle. The above proposed two stage substructural synthesis methodology yields the position-dependent dynamic response at the TCP. A face-mill cutter of 50 mm diameter with an overhang of 70 mm from the spindle nose with a CAT40 type tool-holder is modeled and its response is obtained from Eqs. (5.3) ? (5.5). Response at the spindle nose obtained from solutions to Eq. (5.1) and Eq. (5.6) is synthesized with the tool-tool-holder response, which leads to the response at the free-end of the tool using Eq. (5.7). 5.3. Model Verification and Validation Model checks are made for representative position of the tool being near to the table, i.e. when the spindle housing is at the bottom of its Z axis stroke, at a distance of -0.4 m from the configuration shown in Figure. 5.2. As a first step, response at the TCP predicted with the substructurally synthesized reduced order model is compared with the full order model results for the case of both models having all joints modeled as rigid. Once verified, joints are modeled and identified, and the updated verified model predicted response is validated against measurements. 5.3.1. Contrasting the Substructurally Synthesized Reduced Model with Full Model All connections are initially modeled as rigid, i.e. with in Eq. (5.7) and for displacement compatibility at the interfaces, with the influence of joints being treated in the next section. A uniform damping ratio of ? = 0.02 has been assumed in simulating all the FRFs. These damping ratios are updated by correlating the model predicting response with measured response, in the next Section. Reduced and full model TCP response is compared in Figure 5.3. In addition to comparing the response at the TCP in Figure 5.3, the normalized relative frequency 82 difference (NRFD) and the modal assurance criterion (MAC) between the reduced and full order models are also compared in Figure 5.4. Response comparisons are made up until 350 Hz, i.e. up to the frequency range of interest. Since the machine is envisaged for heavy-duty milling of steels - which are cut at lower spindle speeds (spindle speed being limited by the recommended range of cutting speeds for steel), only the lower frequency structural response is of particular interest to us. The stability of heavy-duty milling is more influenced by the lower frequency structural modes of the larger parts of the machine than the higher frequency tool-tool holder modes, which are in turn damped out by the cutting process in the lower speed zones. Figure 5.3 Full order and reduced order model TCP FRFs comparisons in X (top) and Y (bottom) directions. As is evident from comparisons in Figure 5.3, for the case of rigid connections, the reduced model quite reasonably approximates the full model response. The slight underestimation of natural frequencies and dynamic stiffness by the substructurally 83 synthesized reduced machine model are thought to be due to the nature of the interpolation constraint formulation during syntheses. The interpolation constraint formulation underestimates contact stiffness at the interface, leading to underestimation of stiffness of the overall assembled machine [111]. Figure 5.4 Normalized relative frequency difference (NRFD) and modal assurance criterion (MAC) comparison between the full order and reduced order machine model Since the response in Figure 5.3 is constructed using the modal vectors (mode shapes) of the respective full order and reduced order models, the MAC comparison in Figure 5.4 better explains the discrepancies observed in correctly approximating the dynamic stiffness of the low frequency mode in the machine X direction. Since the full and reduced synthesized 84 models have been modeled differently with different mesh sizes, they do not have the same nodal distribution; and a MAC comparison cannot be made for all DOFs of the reduced model. Hence, MAC comparisons in Figure 5.4 are limited to only a subset of the nodes which are taken as those which can represent the deformation shape of the full machine. These nodal DOFs are chosen as the nodes corresponding to the corners of the guide way top surfaces on the column substructure, the nodes corresponding to the Z-axis feed-drive assembly, and, the spindle nodes. As evident in Figure 5.4, though the diagonal of the MAC for some modes is near unity ? suggesting a good correlation for these modes, there is also a scatter of a few off-diagonal terms which are near unity. This suggests a poor correlation of these modes. The off-diagonal terms observed in Figure 5.4 may be attributed more to the errors due to the synthesis process than due to the reduction process; since as was already shown in Section 3.6.1, the substructural reduced and full models were well correlated with near uniform unity diagonal dominated MAC values. Response comparisons in Figure 5.4 are for the first 15 low-frequency modes within the frequency range of interest, i.e. up until 350 Hz, and not only for the modes observed at the TCP response in Figure 5.3; hence, caution should be exercised in correlating modes in Figure 5.3 with that in Figure 5.4. The errors in frequencies between the reduced and full model fluctuate between -10-15% depending on the mode. A comparison of the natural frequencies in Figure 5.4 also suggests that for the case of a higher frequency error between the reduced and full order model counterpart, the corresponding mode shapes are not well correlated, and have significant off-diagonal terms. Since the objective of the machine tool designer is to maximize the dynamic stiffness, and the reduced model underestimates the dynamic stiffness of some modes, the reduced models may well provide the designer with greater leeway in modifying the design, by erring on the side of caution. Overall however, even though the synthesized reduced model is ~1/25th the size of the full model, yet it gives reasonably close results; and importantly, the reduced model leads to considerable simulation time savings for the designer; taking ~1 minute/ position as compared to ~12 hours/ position for the full model (Intel? i7 2.67 GHz processor with 9 GB RAM) thereby facilitating position-dependent analyses. The verified reduced model is validated against measurements by considering joint characteristics in the next section. 85 5.3.2. Comparing Model Predicted Response and Measured Response Model predicted response is updated by modeling and identifying the joint characteristics as discussed in Chapter 4. The test setup to identify modal damping on the machine being modeled is shown in Figure 5.5. An instrumented impact hammer is used to excite the machine at the spindle nose and a laser vibrometer is used to measure the displacement response. The identified modal parameters (modal damping levels) are subsequently used to update the model predicted response for the lower frequency modes. Joint characteristics for the higher frequency tool-tool-holder modes are identified with the inverse receptance coupling approach as discussed in Chapter 4 and Section 5.2.3. Figure 5.5 Experimental setup to measure TCP FRFs on machine. The predicted response for the updated model obtained by solving Eq. (5.7) with is compared with measured response and with the model in which all joints are treated as rigid in Figure 5.6, for the frequency range of interest, i.e. up to 350 Hz. Measured and model predicted behavior of the dominant modes is also compared in Table 5.2, which also lists the errors. 86 Figure 5.6 Comparison of model predicted response (with/ without joint effects) with measured response As is evident from Figure 5.6 and Table 5.2, the model predicted response, with joint dynamics incorporated, is better able to approximate measured behavior as compared to the rigid model, especially in the machine Y direction. The errors in both models are still considerable, especially for response in the X direction. These errors are carried over from response comparisons at the spindle nose in Chapter 4, Section 4.3.3; and may be attributed to the same reasons outlined before, i.e. due to modeling simplifications and due to the difficulty in approximating contacting interfaces which are idealized as being connected with linear springs. To a lesser extent, the errors may also be attributed to the reduction and synthesis process. The results and trends reported in Table 5.2 are for response at the TCP and are for the reduced model with/without joints, hence are slightly different than those reported in Table 4.2, in which the response was compared at the spindle nose for the full machine models with/without joints. 87 Table 5.2 Comparison of dominant low frequency modes for model predicted response (with/ without joint effects) with measured response at the TCP; where is the natural frequency, is the dynamic stiffness, is the error in natural frequency prediction and is the error in dynamic stiffness prediction Measured Model with joints Rigid model X Mode # [Hz] [N/?m] [Hz] [%] [N/?m] [%] [Hz] [%] [N/?m] [%] 1 36 17 34 -6 9 -90 42 14 6.8 -150 2 97 9.5 60 -60 3.6 -150 75 -30 3.7 -155 3 130 12.5 134 3 22 40 128 -2 25 50 Y 1 26 8 24 -7 6 -33 41 37 14 42 Any error in the prediction of dynamic stiffness at the TCP would lead to errors in chatter stability prediction, and, thence the model would not find much use to establish experimental guidelines for recommending stable cutting parameters. If the model was intended thus, a more detailed and accurate model may be necessary. However, since the purpose of this model is to guide engineering design decisions and to facilitate rapid evaluation and optimization of the machine?s position-dependent dynamic behavior, these errors though high are deemed acceptable at this stage of development. Moreover, since the dynamic stiffness is underestimated, the model errs on the side of caution, i.e. it is more flexible than the physical test structure; leading to potentially optimizing a design concept for a higher stiffness with a higher factor of safety. The model response is simulated and updated with joint characteristics at other positions as well to evaluate the machine?s position-dependent dynamic behavior. 5.4. Position-dependent Dynamic Response at Tool Point TCP FRFs are simulated with the substructurally synthesized reduced order machine model at three different positions of the machine. The top position is the configuration shown in Figure 5.2; the mid and bottom positions are when the tool has moved in the Z-direction by an amount of -0.2 m; and -0.4 m respectively. The stability of the machining process is primarily influenced by the low-frequency structural modes such as column, table, spindle 88 housing, and the spindle shaft; hence X and Y directional direct TCP FRFs are compared up until 350 Hz in Figure 5.7. Higher frequency tool and tool-holder modes are usually not related to the design of the machine tool structure since they are specific to machining application and are more local in nature, i.e. they do not exhibit strong position-dependency. As evident in Figure 5.7, the global modes corresponding to the column (20 to 100 Hz) exhibit stronger position-dependency as compared to the spindle housing modes (100 to 350 Hz). The low-frequency column mode in the X direction at ~34 Hz varies by ~100% in dynamic stiffness, whereas the second X directional column bending mode at ~60 Hz varies by up to 8% in natural frequencies. The dominant low-frequency column mode in the Y direction at ~24 Hz varies by as much as 12% in natural frequencies and ~95% in dynamic stiffness over the full Z stroke of the machine. Figure 5.7 Comparison of reduced model TCP FRFs at three different tool positions: top, mid, and bottom 89 To better understand the low frequency modes and the shape associated with these modes, the first dominant low-frequency mode in each machine direction for the case when the headstock is at the top position are shown in Figure 5.8. The full model shown in Figure 5.8 (a) is represented by only the interface DOFs in Figure 5.8 (b). Mode shapes are shown by overlaying the deformed configuration over the un-deformed configuration. The first mode at ~24 Hz corresponds to a global column bending mode in the YZ plane, in Figure 5.8 (c), while the second dominant mode at ~34 Hz corresponds to the global column bending mode in the XZ plane, in Figure 5.8 (d). These mode shapes do not change as a function of position; position influencing only modal (vibration) amplitudes and frequencies. Figure 5.8 (a) Full order model (b) Reduced model with only interface DOFs, (c) Mode shape of column bending in YZ plane; (d) Mode shape of column bending in XZ plane. 90 Influence of the strong position-dependent dynamics of the verified and validated (within acceptable errors) substructurally synthesized reduced machine model on the performance and productivity levels within the entire work volume will be investigated in the next section. 5.5. Evaluation of Machine Tool Performance 5.5.1. Material Removal Rates and Position-dependent Stability When the structural dynamics of the machine vary within the machine?s work space, the chatter stability and the resulting limits on the material removal rates vary as well. The variation of stability is demonstrated here by considering face milling of AISI 4340 steel with 80% engagement, which is treated as the target application for the envisaged machine. A minimum speed and feed-direction independent stable depth of cut of 3 mm is targeted in the whole working range of the machine. Effects of position-dependent directional compliances on machining stability are investigated by generating feed-direction-dependent absolute machining stability charts. 5.5.2. Oriented FRFs and Feed-direction-dependent Stability A multi-degree-of-freedom machine tool may chatter in any of its dominant modes. The effect of each mode is determined by its dynamic characteristics; whether or not the mode is aligned with the principal machine directions ; and, the direction of feed. Consider a machine tool with horizontal and vertical axes; and let the tool be travelling in the feed direction with an angular orientation as shown in Figure 5.9. Additionally, if the structural mode ? is not aligned with the principal machine tool directions, it would be oriented with an angular distance as also shown in Figure 5.9. 91 Figure 5.9 Orientation of machine mode ( ? ) in the machine tool axes ( ) and the feed axes The stability of the milling system is determined using a modal model of the machine and the following characteristic equation [122]: [ ] [ ] (5.8) where ( ) (5.9) is the eigenvalue of the characteristic equation, and are its real and imaginary parts; is the number of teeth on the cutter; is the cutting force coefficient of the material being cut; is the axial depth of cut; is the chatter frequency; and, is the tooth passing period. The oriented directional matrix within Eq. (16) is expressed as: [ ] [ ] [ ] (5.10) where [ ] is the matrix of the average direction factors, determined as in [122]; and is the 2D transfer function matrix at the TCP in the feed plane, expressed as: [ ] [ ] (5.11) where and are the direct transfer functions in the and feed directions, and where and are the cross transfer functions. In order to obtain the feed-plane transfer function matrix, , the principal modes assumed to be parallel to the directions (i.e. for ) are projected onto the feed ( ) directions with the transformation [103]: 92 { } [ ( )] { } (5.12) where ( ) [ ( ) ( ) ( ) ( )] is the transformation matrix. Similarly, the cutting forces when transformed, become: { } [ ] { } or { } [ ] { }. (5.13) The vibrations in directions: { } [ ] { } { } [ ] { } (5.14) along with Eq. (5.13) when substituted into Eq. (5.12), lead to the vibrations in the feed directions to be expressed as: { } [ ][ ][ ] { } (5.15) wherein [ ][ ][ ] [ ]. Since the principal mode ? is aligned with the machine tool principal axes , which are orthogonal, the cross terms are negligible, i.e. . However, due to the projection of modes in the feed directions, . Representing again the characteristic equation in Eq. (5.8) as a quadratic function in , we get: (5.16) where the coefficients and including the effects of the cross terms from the feed plane FRF matrix are: ( ) ( ); . (5.17) The absolute speed independent limiting depth of cut, described by the parameters in Eq. (5.8) and Eq. (5.10), for different feed directions (0-360?) may be analytically determined as: [ ( ) ] (5.18) 93 Solution to Eq. (5.18) for each tool position within the work volume results in speed independent absolute minimum stable depths of cut which vary with feed directions in proportion to the magnitude of projections of the modes in that direction. 5.5.3. Process Stability and Modes Limiting Productivity Feed-directional stability was simulated for face milling AISI 4340 steel with 5; 3000 MPa; and the radial coefficient, 0.24. Results for three tool positions along with the corresponding chatter frequencies are shown in Figure 5.10. In Figure 5.10 (a), the regions inside the stability curve envelopes are stable, and those outside the process stability threshold, are unstable. The absolute limiting depth of cut is plotted radially, while the machining (feed) directions are plotted circumferentially. These feed-direction dependent absolute stability charts are a departure from the traditional stability lobes [122] in a way that they represent only the absolute minimum stability at a particular position and feed direction and are not concerned with finding the stability lobes which recommend selection of cutting parameters based on spindle speed. The stability curves are symmetric about 45? and 135? respectively, i.e. in directions where the mode under consideration halves the angle between the principal directions. The shape and envelope of the stability boundaries are a strong function of the engagement conditions and position-varying strengths of directional modes. An ideal machine design must have circular stability boundary, i.e. identical depth of cut limits in all directions, which is only possible if the structural dynamics of the machine are uniform in all feed directions. 94 Figure 5.10 (a) Feed directional stability at three different tool positions: top, mid, and bottom, and, (b) Corresponding chatter frequencies as a function of feed direction. Minimum limits of ~1.65 mm at the 50? and 232? feed orientations are observed when the tool is at the bottom position; and, ~1.17 mm at the 52? and 232? feed orientations for when the tool is at the mid position. The absolute minimum limit occurs however when the tool is at the top position of the Z stroke; being ~1 mm at the 52? and 232? feed orientations. This is mainly due to the fact that the dynamic stiffness of modes in the X and Y directions is similar 95 in magnitude and is a minimum at the top position (see Figure 5.7). The feed-directional chatter spectrum in Figure 5.10 (b) also shows that chatter occurs between ~15-65 Hz ? depending on the feed direction, which is near to the low-frequency dominant column bending modes in the XZ and YZ planes (see Figures 5.7-5.8). Having identified the column mode(s) as the cause of limited productivity due to chatter, they need to be stiffened. 5.6. Structural Design Modifications to Meet Target Productivity In order to achieve the target productivity goal, the dynamic stiffness of the modes limiting productivity needs to be increased. The dynamic stiffness being a function of the modal stiffness and the damping coefficient, , an increase may be brought about by one or more of the following strategies: designing structures with a more rigid construction, i.e. by increasing mass and/or stiffness of the substructures; and/or with damping modifications - active or passive. Methods of mass addition/reduction will shift the natural frequencies of vibration, and may result in reduction of vibration amplitude at certain frequencies. However, chatter may still occur at the shifted frequency, since the forcing function may vary as a function of material being cut; workpiece engagement; and cutting speeds. Moreover, mass addition is contrary to the objective of this work, which aims at simultaneous reduction of moving masses. Methods for tuning through passive damping, relying on better energy transfer to the ground through mounting pads etc. though favorable, is outside the scope of the present discussion. Active damping though effective, requires additional external actuators and adds to the complexity of the system, hence is not the preferred method for dynamic stiffness modification in the present discussion. With the forgoing discussion, structural modification by stiffening is preferred in this work; due to its relative ease in implementation at the design stage. The dominant modes being the bending and torsion of the column (see Figures 5.7-5.8), the column substructure is stiffened through addition of flange stiffeners and internal cross ribs. These modifications shown in Figure 5.11 result in increase in mass of the column by ~15%, which is contrary to the objective of development of light-weight machine tools. Hence, a subsequent topology optimization is carried out to optimize material distribution within the column design space while minimizing its volume. Topology optimization has been shown to be a useful tool in 96 supporting design decisions [7,123] is carried out on the Hypermesh? [124] platform in the present investigations. Figure 5.11 (a) Original, Modified, Optimized, and Smoothened Column models; (b) Comparison of X direction TCP FRFs for original and modified machine models; and, (c) Comparison of Y direction TCP FRFs for original and modified machine models. The design space for the column is all the parametric space other than the contacting interfaces. A new design concept is found by subjecting the TCP to loads representing 97 machining forces of 3000 N while respecting constraints on compliance. Optimization investigations have been conducted at all three positions of the spindle housing described earlier, and the final design concept is shown in Figure 5.11 (a). The new design concept for the column has the same mass as that of the original column model, i.e. a ~15% reduction in mass from the stiffened model was achieved through optimization. This design concept is smoothened into a suitable geometry-based CAD model for further analysis. The procedure given in the flowchart in Figure 5.1 is followed again with the new column concept to obtain the new synthesized reduced machine tool model. The reduced model?s dynamic response at the TCP for the original, stiffened and subsequently optimized machine model is compared in Figure 5.11 (b, c) for a representative condition when the tool is at the top position; and, with all connections assumed to be rigid. For comparisons, a uniform damping ratio of ? = 0.02 has been assumed. Joint characteristics are incorporated through modeling and identification, as done before, after these initial comparisons and before simulating position-dependent dynamics and stability for the new machine concept. An increase in dynamic stiffness by ~20% for each of the first bending modes in the XZ and YZ planes; and by ~30% for the mode at ~90 Hz is evident from TCP FRF comparisons in Figure 5.11 (b,c), in which comparisons are limited to 120 Hz, i.e. the frequency range affected by the structural modifications. The increase in the natural frequencies of the low-frequency modes is due to stiffening. The dynamic response with the optimized topology is equivalent to the response with the stiffened structure with the advantage of having less overall structural mass. The static stiffness for the modified machine design now ranges from 20-22 N/?m (depending on tool position), an increase of ~40% over the previous range of 14-16 N/?m. 5.6.1. Improved Machine Tool Performance Stability was simulated with the same tool and cutting conditions as before, but with the modified machine response; and the stability results are compared with the values of the original machine model for the three positions as before, see Figure 5.12. 98 Figure 5.12 Feed-directional stability comparison for original and modified (optimized) machine model at two different tool positions: (a) bottom, (b) mid, (c) top As is evident from Figure 5.12, the new absolute stability limit at the top position is ~1.5 mm at the 52? and 232? feed-orientations, an increase of ~50% over the limit achievable with the original machine model. At the mid position, the new absolute limit is ~1.7 mm, an improvement by ~45% over the earlier ~1.17 mm limit; and at the bottom position, the new absolute limit is ~2.1 mm, an increase of ~30% over the earlier ~1.65 mm limit. The shape and envelope of the stability boundaries are different than the original design due to the change in position-varying directional compliances for the modified structure. Structural 99 modifications result in an increase in the stability envelope leading to higher absolute stability limit for a wider range of feed directions, even exceeding the target productivity levels of 3 mm at certain feed directions for all tool positions; which is a significant improvement over the earlier design. Further improvements, if necessary, to meet the target 3 mm depth of cut at all working positions and feed-orientations are possible by modifying (stiffening) other substructures, or, by integrating active/ passive damping mechanisms ? all of which are facilitated by the rapid iterative procedure developed in this Chapter. 5.7. Summary This Chapter demonstrates the feasibility of a virtual engineering approach to rapidly evaluate and optimize the performance of machine tools through design modifications to improve their dynamic stiffness. By characterizing the position-dependent process-machine interactions, the machine with improved dynamic stiffness is shown to better achieve targeted productivity. The position-dependent substructurally synthesized reduced order machine model efficiently and accurately simulates dynamic response; taking ~1 minute/ position as compared to ~12 hours/ position for the full model. The model facilitates rapid investigation of design alternatives compared to the time consuming full FE models used presently. For a defined set of representative milling operations and target productivity levels for which the machine is being envisaged, the effects of position and feed-direction-dependent compliances and machining stability are investigated. It is shown that the proposed reduced order substructural model allows rapid redesign and analysis of components which limit the productivity due to low compliance. Though the model predicted response underestimates the dynamic stiffness as compared to measured behavior due to difficulties in accurately modeling the joint characteristics, it errs on the side of caution; i.e. the resultant stability maps generated with model predicted response envelope smaller areas and have a lesser absolute minimum stability limit than if the stability maps were to be generated with measured response. This translates to a higher factor of safety in the design process, resulting in a potentially stiffer design than would be necessary for achieving the target dynamic stiffness and productivity levels. 100 Even though the design methodology was applied to a specific kinematic configuration and for a given set of machining operations, the methodology developed is generic and may be extended to virtual evaluation of any machine tool ? as is done in the next chapter, in which the position-dependent dynamics and stability of a serial-parallel kinematic machine tool are assessed. The proposed virtual machine concept facilitates total system simulation, including the dynamics of the control and mechanical parts of the machine; making it possible to assess design suggestions and to perform system optimization at an early stage of development. 101 Chapter 6 Position-dependent Dynamics and Stability Analyses of Serial-Parallel Kinematic Machines1 6.1. Overview This Chapter extends the substructurally synthesized reduced machine modeling scheme developed in Chapter 3 to model and evaluate the position-dependent dynamics for a hybrid serial-parallel kinematic machine tool developed for machining large workpieces in the die and mold industry. Such machines must necessarily envelope large work volumes and deliver consistently high dynamic stiffness over the entire work volume. This is difficult to achieve with traditional serial machine tool construction, since higher stiffness necessarily involves bigger structures, which require bigger, more powerful and costly feed drive motors. Alternatively, nimbler parallel kinematic machine tools with higher dynamic stiffness capabilities have been developed; however, these suffer from strong position-dependent dynamics, leading to sluggish performance over large work volumes. To deliver improved dynamic performance over the whole working range of the machine, a prototype hierarchical hybrid serial-parallel scissor kinematic machine has been developed at the Fraunhofer IWU [125], see Figure 6.1. Motion is split between the serial parts that envelope large work volumes and the scissor kinematic arrangement that allows delivery of higher dynamic stiffness at the tool center point (TCP) with higher acceleration capabilities [126]. Though the scissor arrangement with struts of fixed lengths is superior to other parallel kinematic arrangements with struts of varying lengths, the scissor kinematics still exhibits strong position-dependent dynamic behavior [51]. The position-varying dynamics lead to position-varying machining stability of the system, which limits the achievable productivity and performance in the whole working range of the machine [57]. The evaluation of the changing stability helps to plan stable machining trajectories by splitting the motion between the serial and parallel parts. In order to exploit the potentially superior dynamic behavior of such a hybrid arrangement it is necessary to consider the position-dependency of the system at the design and development stage. 1 Portions of this Chapter have been published as Law, M., Ihlenfeldt, S., Wabner, S., Altintas, Y. and Neugebauer, R., 2013, Position-dependent dynamics and stability of serial-parallel kinematic machines, CIRP Annals - Manufacturing Technology, Vol. 62 [4]. 102 Figure 6.1. Hybrid serial-parallel scissor kinematic machine. To facilitate rapid assessment of position-dependency for such hybrid machines, the modeling processes developed in Chapter 3 are extended here, leading to a novel method to model the parallel struts undergoing rotation in Section 6.2. The machine model is verified and validated against measurements in Section 6.3, followed by simulating and subsequently validating the machine?s position-dependent behavior in Section 6.4. The reduced multibody machine model is used to simulate position-dependent stability maps in Section 6.5 so as to assess stable material removal rates over the entire work volume. 6.2. Position-dependent Multibody Dynamic Machine Model The machine is treated as an assembly of the major flexible substructures as shown in Figure 6.1. Motion in the XY plane is achieved either by four linear motors driving the Y-slides, or, through actuating the scissor-kinematic arrangement as shown in Figure 6.2. Struts at either side of the platform are driven by action of four independent linear motors, whose planar motion results in rotation of the struts about their pivot points. Independent control of the linear motors may inadvertently cause rotation of the tool platform about the Z axis. To avoid this, there are three struts on the one side and two on the other [125]. Tool motion in 103 the XZ plane is achieved by two ball-screw feed drives. For details on the kinematics of the machine, see [125,126]. Figure 6.2. Kinematic scheme for the scissor kinematics. During positioning, the substructures undergo simultaneous translation and rotation, which leads to position dependent structural dynamics. First, each major substructure is modeled independently in the FE environment; and is exported, after convergence checks to the MATLAB environment for reduction and synthesis with the parallel struts. The parallel struts in turn are modeled with Timoshenko beam elements which are oriented to the appropriate configuration using transformation matrices prior to assembly with the main structural components. The details of general mathematical model based on dynamic substructuring remain the same as discussed earlier in Chapters 3, 4, and 5; and, only the principles that are specific to the hybrid serial-parallel kinematic machine are given here. 104 6.2.1. Substructure Model Reduction FE models of the each of the major structural components shown in Figure 6.1, i.e. the two Y-slides, four guide-blocks, the Z-carriage, the Z-slide and the two ball-screw feed-drive models for the Z axis, are all modeled as discussed in Chapter 3, Section 3.2. The structural components were assigned material properties of structural steel as: modulus of Elasticity of 210 GPa; density of 7800 kg/m3; and, Poisson?s ratio of 0.3. The two Y-slides are made of a special sandwich foam in which aluminum foam is sandwiched between two steel plates [127]. Modeling such detail within the FE environment is difficult and cumbersome. Moreover, a detailed sandwich model within FE results in a very large order model which poses difficulties in handling within the MATLAB environment during the reduction process. Hence, material properties for these Y-slides are selected such that the simplified model response characteristics match the detailed model characteristics. Each of the substructural models is exported after convergence checks to the MATLAB environment for reduction and synthesis as discussed in Chapter 3, Sections 3.3-3.4. These structural models are synthesized with the parallel struts which are modeled with Timoshenko beam elements as discussed in the next section. 6.2.2. Modeling Parallel Kinematics Each of the five parallel struts is of a hollow box type construction and may also be modeled as the other structural components, i.e. with solid tetrahedral 10-noded elements. The reduced model substructural synthesis approach is based on the synthesis of the individual substructures that have position-independent response to obtain position-dependent response. This is achieved by first aligning each of the individual substructures to the correct position and subsequently enforcing displacement compatibility at the substructural interfaces. This is rather straightforward for a machine that has substructures undergoing translational motion. However, in the present situation, the parallel struts undergo large rotation about the pivot points due to the translating guide-blocks. If the strut was to be modeled with solid elements, it would not be possible to rotate the meshed strut to the proper orientation outside of the FE environment before synthesis. Hence, an alternate approach is proposed to model the struts, i.e. one in which each of the parallel strut is modeled with Timoshenko beam elements. Beam elements have the distinct advantage that they may be 105 rotated rather straightforwardly to obtain the desired orientation with the help of a rotation (transformation) matrix. Each beam element has six DOFs on each of its two nodes; three translational ( ) and three rotations ( ) ? as shown in Figure 6.3. A two stage transformation is carried out at the elemental level: i.e. first a rotation about the X axis ( ) to bring the strut to the correct inclination as in Figure 6.1; and, a subsequent rotation about the Z axis ( ) to make the correct orientation as in Figure 6.2. Transformations are based on kinematic relationships. Figure 6.3. Oriented parallel strut modelled with Timoshenko beams For a given tool position, the platform subtends different angles at either side ( ), except for the tool position for when the platform is at the center of the X2 axis, in which case, . The inclination ( ) for each of the struts is fixed, i.e. it does not change as a function of tool position. The formulation discussed here is for the general case for any of the five parallel struts. 106 The transformed structural matrices { } for each of the parallel struts at the elemental level are expressed as: (6.1) where { } are the elemental matrices in YZ plane; and is the two stage transformation matrix expressed as: (6.2) where [ ] [ ] (6.3) wherein the nodal operator at each node j is expressed as: [ ] [ ] (6.4) Once transformed to the appropriate inclination and orientation, each strut is assembled with a total of elements as shown in Figure 6.3. In this manner each of the parallel struts may be modeled independently, followed by bringing it to the correct orientation and inclination based on the transformation matrices (Eqs. (6.1) ? (6.4)) prior to synthesis with the other structural components. 6.2.2.1 Response Comparisons for Parallel Strut Modeled with Solid and Beam Elements As a first check, to verify if indeed the strut may be approximated with beam elements, response comparisons are made for the case of the parallel strut modeled with solid and beam elements in the free-free configuration of the strut. The parallel strut modeled with solid 107 elements is shown in Figure 6.4 and the free-free FRFs between the connection ends are compared in Figure 6.5. Figure 6.4 Solid Model of the parallel strut Figure 6.5 Comparison of free-free response at free (connection) ends of the parallel strut modeled with different element types A uniform damping ratio of ? = 0.02 has been assumed in simulating the FRFs. As is evident from comparisons in Figure 6.5, the response with the beam is able to reasonably capture the first two dominant modes within the frequency range of interest. 108 6.2.3. Reduced Model Substructural Synthesis Having verified that the parallel strut may be represented by beam elements, the assembled structural matrices of each of the parallel struts are synthesized at either end with the reduced substructural models. The undamped equation of motion of the synthesized reduced machine model at a particular position is: [ ] { ? ? ? ? ? ? } [ ] { } { } (6.5) where { } are the reduced substructural mass and stiffness matrices obtained with the improved reduced order modeling scheme as discussed in Chapter 3; and, { } are the oriented assembled matrices corresponding to each of the five struts. The subscripts and correspond to the Z-slide (including the tool, tool-holder and spindle); the Z-carriage; each of the five parallel struts; the four guide-blocks; and, both Y-slides respectively. ( ) is a displacement operator whose coefficients are obtained by ensuring displacement compatibility at the contacting interfaces with an interpolation multipoint constraint equation formulation, as in Chapter 3. in Eq. (6.5) represent a discrete set of Lagrange multipliers which explicitly satisfy the constraint formulation. The synthesized model enables prediction of dynamic response at the TCP as one component changes its position/orientation relative to another by solving the eigenvalue problem form of Eq. (6.5). The overall sequence of operations to obtain the position-dependent response is as shown schematically in Figure 6.6. Each new tool position results in new oriented structural matrices for the parallel struts in Eq. (6.5). Since the rotation of the struts causes simultaneously changing boundary conditions on the sliding interfaces, the 109 displacement operator in Eq. (6.5) gets suitably updated, by instantaneously coupling/de-coupling the corresponding nodes on the interfaces. Figure 6.6 Overall flow chart showing the sequence of operations to obtain substructurally synthesized position-dependent response for the hybrid machine 6.3. Model Verification and Validation Model validation requires joints to be modeled in the FE environment; which is non-trivial due to the dependence of joint characteristics on parameters like contact surface conditions, friction and damping. To simplify joint modeling, the same modeling strategy as adopted in Chapter 4 and 5 is followed here, i.e. the contacting interfaces are idealized as being connected by spring elements for which the equivalent contact stiffness values are obtained from manufacturers? catalogues for the rolling interfaces; and, are taken as the corresponding bolt?s stiffness for the bolted interfaces. Damping at the joints is incorporated by updating model predicted response with measured modal damping ratios. Updated model predicted (full, reduced) frequency response functions (FRFs) are compared with measured FRFs in Figure 6.7; in which the dominant modes are labeled. FRFs 110 at the Z-slide surface are compared at the representative position of the tool at the mid position of the X2 axis stroke and bottom position of the Z axis stroke. FRFs are compared up until only the frequency range of interest, i.e. up to 60 Hz. For a large machine tool as this, only the lower frequency structural response is of particular interest to us since the stability of heavy-duty milling (typical operation envisaged for this machine) is more influenced by the lower frequency structural modes that belong to the larger parts of the machine than the higher frequency tool-tool holder modes. Figure 6.7 FRF comparisons between measured, reduced and full model. As evident in Figure 6.7, model predicted response reasonably captures the trend of the measured response characteristics. Moreover, a good match between the number of measured and simulated modes is also observed. The substructurally synthesized reduced machine 111 models underestimate the natural frequencies as well as the dynamic stiffness of the dominant modes. This is for the same reasons as before, i.e. due to the nature of the substructural synthesis process in which the approximate models of surface interactions at the contacting interfaces underestimates the contact stiffness, thereby underestimating both the natural frequencies and the dynamic stiffness. Errors may also be attributed to the difficulties in modeling the joints using simplified ?whole-joint? approximations, and partially due to modeling simplifications. As highlighted earlier in Section 5.3.2 in Chapter 5, if the model was intended to provide experimental guidelines for recommending stable cutting parameters, the model may find limited applicability, and a more detailed and accurate model may be necessary. However, since the purpose of this model is to rapidly evaluate the potential superiority of a hybrid serial-parallel kinematic in terms of its changing dynamic behavior, these errors are deemed acceptable at this stage of development. Moreover, since the reduced model size at ~13,500 DOFs is ~1/30th the size of the full model, which has ~400,000 DOFs; it leads to considerable simulation time savings for the designer; taking ~4 minutes/ position as compared to ~20 hours/ position for the full model (Intel? i3-380M processor with 4 GB RAM) thereby facilitating position-dependent analyses. 6.4. Position-dependent Dynamics Influence of the orientation of the parallel struts and of the Z position of the platform on the dynamic behavior at the TCP, i.e. on the natural frequencies ( ) and dynamic stiffness ( ) of the low frequency dominant modes in the X and Y directions is evident from comparisons in Figure 6.8; in which position-dependent response is predicted over the whole working range of the X2 axis (? 150 mm) and the Z axis stroke (? 900 mm). Weak position-dependency was observed over the Y stroke of the machine. Measured position-dependent response is overlaid over predicted response. The trend in errors observed at a single position in Figure 6.7 is carried over to the different TCP positions within the work volume, and as evident in Figure 6.8, though there exist some discrepancies, model predictions reasonably approximate measured behavior. The response characteristics for the dominant modes in machine X and Y directions, i.e. the stiffness as well as the natural frequencies are observed to change more strongly over the 112 Z-stroke of the machine than over the X2 axis stroke. Interestingly, the response characteristics are not symmetric about the neutral position ( ) of the tool platform along the X2 axis. This may be explained by the fact that there are three struts on one side of the platform, and two on the other, making one side potentially stiffer than the other. Model predicted response is able to capture this behavior. Moreover, the kinematic design results in response in the machine Y direction to vary more than in the X direction ? which is also captured by the model. Figure 6.8 Measured (?) and predicted (surface plots) position-dependent behaviour of the machine in the X2Z plane. 113 6.5. Evaluation of Position and Feed-direction-dependent Stability The variation of stability due to position-varying structural dynamics is demonstrated here by considering face milling of AISI 4340 steel with 80% engagement, which is again treated as the target application for the envisaged machine. As observed earlier in Chapter 5, wherein the effects of position-dependent directional compliances on machining stability were investigated by generating feed-direction-dependent absolute machining stability charts; the same procedure is adopted here in determining the absolute speed independent limiting depth of cut for different feed directions (0-360?) using the following: [ ( ) ] (6.6) where all symbols are already defined in Section 5.5.2, Chapter 5. Solution to Eq. (6.6) for each tool position within the work volume results in speed independent absolute minimum stable depths of cut which vary across feed directions in proportion to the magnitude of projections of the modes in that direction. The tool and work piece parameters considered are: face mill of diameter 100 mm with 80 mm overhang from the spindle nose with 10; 3000 MPa; and the radial coefficient, 0.24. Absolute feed-direction-dependent stability limits were generated for all platform positions in Figure 6.8; and, the resulting stability envelope covering the whole work volume is shown in Figure 6.9. The region inside the stability envelope is stable, and that outside, is unstable. The absolute limiting depth of cut is plotted radially, while the machining (feed) directions are plotted circumferentially. The stability threshold in Figure 6.9 for every Z-position of the platform at each feed direction is defined as the minimum of the stable depths of cut of all tool positions within the X2Y plane. The shape and envelope of the stability boundaries are a strong function of the engagement conditions and position-varying strengths of directional modes. Minimum limits of ~0.5 mm at the 52? and 232? feed orientations are observed when the tool is at the bottom position; and, ~0.9 mm at the 52? and 232? feed orientations when the tool is at the top position. 114 Figure 6.9 Position-varying feed-direction-dependent stability over the whole work volume. The stability map in Figure 6.9 can be used to guide design modifications to improve material removal rates as was done in Chapter 5; and/ or can also serve as guidelines for planning of dynamically stable machining trajectories. 6.6. Summary The influence of position-varying dynamics of a serial-parallel scissor kinematic machine on the dynamic stability of a milling process is investigated in this Chapter. Position-dependency is modeled based on syntheses of reduced substructures which have position-invariant response. A novel method to model the parallel struts undergoing simultaneous translation and rotation is proposed based on representing the struts using Timoshenko beams. The struts are rotated and oriented to the desired configurations with the aid of rotational operators prior to synthesis with the other reduced substructures, thus facilitating 115 efficient modeling and prediction of position-dependent response without having to use the presently used time consuming full order FE models. The model is verified against its full model counterpart and is validated against measurements. Dynamic behavior is shown to be strongly dependent on the orientation of the parallel struts and on the Z-position of the platform in the work-volume. The influence of changing dynamics on the machining process is assessed by constructing stability maps which establish the absolute minimum stable depths of cuts for any given position of the tool in the active work volume. These stability maps serve as guidelines to evaluate the potentially superior dynamic behavior of such hybrid serial-parallel kinematic machine over other conventional serial and/ or parallel kinematic machine tools. 116 Chapter 7 Conclusions and Future Research Directions 7.1. Conclusions Machine tools exhibit position-dependent dynamic behavior as the cutting tool travels along the tool path, which changes the stability and the productive cutting conditions within the machine work volume. A generalized computationally efficient substructurally synthesized reduced order model is developed that facilitates rapid evaluation and optimization of the machine tool?s position-dependent dynamic performance without having to use traditionally used time consuming full FE models. The generalized formulation was used to evaluate the position-dependent behavior of two separate machine tools ? one with a serial kinematic configuration, and the other with hybrid serial-parallel kinematics. The bottom-up approach allows substructures that are modeled independently and that have position-invariant response to be synthesized at the desired configuration to yield position-dependent response. Improved reduced models that are developed are shown to be more accurate than the commonly used standard CMS methods. The reduced substructural models are synthesized at their contacting interfaces using novel adaptations of constraint formulations that tolerate mesh incompatibility at the contacting interfaces during synthesis. This allows for a generalized, simple and efficient synthesis procedure that facilitates modularity in the design and synthesis process. The substructurally synthesized reduced machine model is verified against full order models. Verified models are subsequently validated against measurement by updating models with joint characteristics that are modeled and identified using a two-stage dynamic substructuring approach. Verified and validated models are subsequently used to demonstrate the feasibility of a virtual engineering approach for designing machine tools to improve their dynamic performance. For a defined set of representative milling operations and target productivity levels, the effects of position and feed-direction-dependent compliances and machining stability are investigated. The influence of changing dynamics on the machining process was assessed by constructing stability maps which establish the absolute minimum stable depths of cuts for any given position of the tool in the active work volume. It is shown that the proposed reduced order substructural model allows rapid redesign and analysis of components which limit the productivity due to low compliance. 117 In summary, main contributions of this thesis are: ? An improved variant of the component mode synthesis method has been developed for model reduction. The proposed model better represents higher order dynamics of the substructure while keeping the order of the reduced model to a minimum by spanning a much wider frequency range with fewer component modes than would be required with standard CMS methods. ? Novel adaptations of constraint formulations are developed that tolerate mesh incompatibility during reduced model substructural synthesis. The constraint formulations describe displacement compatibility between adjacent substructures with sets of algebraic equations, which are updated to account for relative motion; thereby efficiently and accurately simulating the position-dependent response. ? A two-stage substructural synthesis procedure to model and identify joint characteristics between various machine tool substructural interfaces is proposed. The response at the spindle nose obtained from the substructurally synthesized reduced machine model is combined with the tool-tool-holder response which in turn is obtained by synthesizing the response of the tool with the tool-holder response using an inverse receptance coupling approach. This two-stage procedure reduces the need to perform several measurements at the spindle nose that was required in earlier studies to obtain translational and rotational FRFs, thereby facilitating rapid investigations of the effect of joint characteristics on machine response. ? A novel position and feed-direction-dependent-process-stability performance criterion is proposed to evaluate the productivity of machine tools in its entire work volume. This criterion also identifies parameters limiting the target productivity levels which are modified (optimized) to meet design targets. ? Influence of the position-varying dynamics of a serial-parallel scissor kinematic machine on the dynamic stability of a milling process was investigated by using a novel method to model the parallel struts undergoing simultaneous translation and rotation based on representing the parallel struts using Timoshenko beams. This thesis presents a generalized and computationally efficient flexible multibody dynamic model that facilitates rapid evaluation and optimization of the position-dependent performance of machine tools. 118 7.2. Future Research Directions The accuracy of the substructural synthesis formulation is based on linear weighted constraint relationships which may suffer from numerical issues due to being dependent on weighting functions. This may be addressed by exploring alternate advanced models of surface interactions based on higher order constraint formulations that are of the form of weighted polynomials and/or are based on a non-uniform rational basis spline (NURBS) representation of the surfaces in contact. Accurate prediction of machine tool structural dynamics depends on the availability of high fidelity joint models that are truly predictive, and which rely less on updating models from experimentally obtained modal damping. To address this, standard libraries and databases could be developed for all the different joint types used in machine tools along with their experimentally obtained damping coefficients. Moreover, for all other things being constant, since the absolute stable depth of cut is directly proportional to the level of modal stiffness and damping and any bounds placed on them, uncertainty analysis for joint modeling/identification may also be carried out to investigate the effects of these uncertainties on chatter stability limits. Based on the position-dependent dynamic models developed in this thesis, control-structure interactive effects could be investigated such as to ensure accurate tracking and positioning of the tool. Moreover, the position-varying plant models could be used in virtual investigations in the development of advanced robust control strategies. This thesis presented a generalized position-dependent multibody dynamic model capable of modeling linear (translational) motion of the tool in the work volume. However, since most high-performance machine tools have a combination of linear and rotational axes ? the model could be extended to also consider the effects of orientation-dependent dynamics on tool point compliance. Though structural level modification followed by topology optimization was preferred in this thesis to improve the machine tool?s dynamic performance, other methods such as passive and /or active means of vibration suppression such as to increase the chatter-free productivity levels of machine tools may also be pursued. 119 Additionally, in-process adjustment of the tool path based on an FE calculated dynamic deformation of the machine may also be explored as an application for this reduced order position-dependent model. The ultimate aim would be to develop a complete virtual mechatronic model of machine tools which would include the position-varying structural dynamics combined with models for process-machine interactions and control-structure interactions. This synthesized mechatronic machine tool system will allow evaluation and optimization of the dynamic performance of the machine tool at the design stage, reducing the need for costly physical prototyping and in turn lead to development of light-weight high productivity eco-efficient machine tools. 120 Bibliography [1] Law M., Phani A. S., and Altintas Y., 2013, ?Position-Dependent Multibody Dynamic Modeling of Machine Tools Based on Improved Reduced Models,? ASME Journal of Manufacturing Science and Engineering, 135(2). [2] Mehrpouya, M, Law, M, Park, C, Park, S, Altintas Y., 2013, ?Joint Dynamics Identification on a Vertical CNC Machine,? 2nd International Conference on Virtual Machining Process Technology, Hamilton, Canada. [3] Law M., Altintas Y., and Phani A. S., 2013, ?Rapid evalution and optimization of machine tools with position-dependent stability,? International Journal of Machine Tools and Manufacture, 68, pp. 81?90. [4] Law M., Ihlenfeldt S., Wabner M., Altintas Y., and Neugebauer R., 2013, ?Position-dependent dynamics and stability of serial-parallel kinematic machines,? CIRP Annals - Manufacturing Technology, 62(1). [5] Erdel B., 2003, High-speed machining, Society of Manufacturing Engineers, Dearborn, MI. [6] Homann B. S., and Thornton A. C., 1998, ?Precision machine design assistant?: A constraint based tool for the design and evaluation of precision machine tool concepts,? Artificial Intelligence for Engineeing Design, Analysis and Manufacturing. [7] Altintas Y., Brecher C., Weck M., and Witt S., 2005, ?Virtual Machine Tool,? CIRP Annals - Manufacturing Technology, 54(2), pp. 115?138. [8] Slocum A., 1992, Precision Machine Design, Prentice-Hall. [9] Zatarian, M, Lejardi, E, Egana F., 1998, ?Modular Synthesis of Machine Tools,? CIRP, 47, pp. 333?336. [10] Reinhart G., and Wei?enberger M., 1999, ? Multibody Simulation of Machine Tools as Mechatronic Systems for Optimization of Motion Dynamics in the Design Process,? Proc. of the IEEE/ASME Int. Conf. on Advanced Intelligent Mechatronics, Atlanta, USA, pp. 605?610. [11] Craig Jr R. R., and Bampton M. C. C., 1968, ?Coupling of substructures for dynamics analyses,? AIAA Journal, 6, pp. 1313?1319. 121 [12] Zhang G. P., Huang Y. M., Shi W. H., and Fu W. P., 2003, ?Predicting dynamic behaviours of a whole machine tool structure based on computer-aided engineering,? International Journal of Machine Tools and Manufacture, 43(7), pp. 699?706. [13] Mottershead J. E., and Friswell M. I., 1993, ?Model updating in structural dynamics: a survey,? Journal of Sound and Vibration, 167, pp. 347?375. [14] Segalman D., Bergman L., and Ewins D. J., 2006, Report on the SNL/NSF International Workshop on Joint Machanics, Arlington, Virginia. [15] Friswell, M I, Mottershead J. E., 1995, Finite Element Model Updating in Structural Dynamics, Kluwer Academic Publishers. [16] Zulaika J. J., Campa F. J., and Lopez de Lacalle L. N., 2011, ?An integrated process?machine approach for designing productive and lightweight milling machines,? International Journal of Machine Tools and Manufacture, 51, pp. 591?604. [17] Reddy, C, Rao S., 1978, ?Automated Optimum Design of Machine Tool Structures for Static Rigidity, Natural Frequencies and Regenerative Chatter Stability,? Journal of Engineering for Industry, 100, pp. 137?146. [18] Rao, S, Grandhi R., 1983, ?Optimum Design of Radial Drilling Machine Structure to Satisfy Static Rigidity and Natural Frequency Requirements,? Journal of Mechanisms, Transmissions, and Automation in Design, 105, pp. 236?241. [19] Yoshimura M., Takeuchi Y., and Hitomi K., 1984, ?Design Optimization of Machine- Tool Structures Considering Manufacturing Cost , Accuracy , and Productivity,? 106(84), pp. 531?537. [20] Koenigsberger F., and Tlusty J., 1970, Machine Tool Structures, Volume 1, Pergamon Press. [21] Catania G., and Mancinelli N., 2011, ?Theoritical-experimental modeling of milling machines for the prediction of chatter vibration,? International Journal of Machine Tools and Manufacture, 51, pp. 339?348. [22] Rahman M., Mansur M. A., and Chua K. H., 1988, ?Evaluation of Advanced Cementitious Composites for Machine-Tool Structures,? CIRP Annals, 37(1), pp. 373?376. 122 [23] Bianchi G., Paolucci F., Van den Braembussche P., Van Brussel H., and Jovane F., 1996, ?Towards virtual engineering in machine tool design,? Annals of CIRP, 45, pp. 381?384. [24] Zaeh, M F, Oertli T., 2004, ?Finite Element Modeling of Ball Screw Feed Drive Systems,? CIRP Annals - Manufacturing Technology, 53(1). [25] Okwudire C. E., and Altintas Y., 2009, ?Hybrid Modeling of Ball Screw Drives With Coupled Axial, Torsional, and Lateral Dynamics,? Journal of Mechanical Design, 131(7), p. 071002. [26] Kamalzadeh A., 2008, ?Precision Control of High Speed Ball Screw Drives,? Ph.D. Thesis, University of Waterloo. [27] Cao Y., and Altintas Y., 2004, ?A general method for the modeling of spindle-bearing systems,? ASME Journal of Mechanical Design, 126, pp. 1089?1104. [28] Kolar P., Sulitka M., and Janota M., 2011, ?Simulation of dynamic properties of a spindle and tool coupled with a machine tool frame,? The International Journal of Advanced Manufacturing Technology, 54, pp. 11?20. [29] Yu, Z, Nakamoto, K, Takeuchi Y., 2011, ?Development of an Interactive Assistance System for Machine Tool Structure Design Considering of Sliding Joint Damping,? International Journal of Automation Technology, 5(5), pp. 722?728. [30] Maj, R, Modica, F, Bianchi G., 2005, ?Machine tools mechatronic analysis,? Proc. IMechE Part B: J. Engineering Manufacture, 220, pp. 345?353. [31] Mori M., Piner Z., Ding K., and Hansel A., 2008, ?Virtual evaluation for machine tool design,? Proceedings of the 9th Biennial ASME Conference on Engineering Systems Design and Analysis, Haifa, Israel. [32] Hung J.-P., Lai Y.-L., Lin C.-Y., and Lo T.-L., 2011, ?Modeling the machining stability of a vertical milling machine under the influence of the preloaded linear guide,? International Journal of Machine Tools and Manufacture, 51(9), pp. 731?739. [33] Jonsson A., Wall J., and Broman G., 2005, ?A virtual machine concept for real-time simulation of machine tool dynamics,? International Journal of Machine Tools and Manufacture, 45, pp. 795?801. 123 [34] Van Brussel H., Sas P., Nemeth I., Fonseca P. D., and Van den Braembussche P., 2001, ?Towards a Mechatronic Compiler,? IEEE/ASME Transactions on Mechatronics, 6(1), pp. 90?105. [35] Neugebauer R., Scheffler C., and Wabner M., 2011, ?Implementation of control elements in FEM calculations of machine tools,? CIRP Journal of Manufacturing Science and Technology, 4(1), pp. 71?79. [36] Zaeh M. F., and Hennauer M., 2011, ?Prediction of the dynamic behavior of machine tools during the design process using mechatronic simulation models based on finite element analysis,? Production Engineering Research and Development, 5, pp. 315?320. [37] Zaeh M. F., and Siedl D., 2007, ?A New Method for Simulation of Machining Performance by Integrating Finite Element and Multi-body Simulation for Machine Tools,? Annals of CIRP, 56, pp. 383?386. [38] Aurich J. C., Biermann D., Blum H., Brecher C., Carstensen C., Denkena B., Klocke F., Kroger M., Steinmann P., and Weinert K., 2009, ?Modelling and simulation of process: machine interaction in grinding,? Production Engineering Research and Development, 3, pp. 111?120. [39] Light, T, Gorlach, I, Wiens G., 2011, ?Dynamic Modelling of a Special Purpose CNC Machine,? Proceedings of the 15th WSEAS international conference on Systems, Wisconsin, pp. 207?213. [40] Okwudire C., 2009, ?Modeling and Control of High Speed Machine Tool Feed Drives,? Ph.D. Thesis, University of British Columbia. [41] Da Silva M. M., Bruls O., Swevers J., Desmet W., and Van Brussel H., 2009, ?Computer-aided integrated design for machines with varying dynamics, ,? Mechanism and Machine Theory , 44, pp. 1733?1745. [42] Sekler P., Voss M., and Verl A., 2012, ?Model-based calculation of the system behavior of machine structures on the control device for vibration avoidance,? The International Journal of Advanced Manufacturing Technology, 58, pp. 1087?1095. [43] Hoher S., and R?ck S., 2011, ?A contribution to the real-time simulation of coupled finite element models of machine tools ? A numerical comparison,? Simulation Modelling Practice and Theory, 19(7), pp. 1627?1639. 124 [44] Cho, M, Kim W., 2002, ?A Coupled Finite Element Analysis of Independently Modeled Substructures by Penalty Frame Method,? KSME International Journal, 16(10), pp. 1201?1210. [45] Aminpour, M A, Ransom, J B, McCmeary S. L., 1995, ?A coupled analysis method for structures with independently modelled finite element subdomains,? International Journal for Numerical Methods in Engineering, 38(21), pp. 3695?3718. [46] Farhat C. G., and Geradin M., 1994, ?On a component mode synthesis method and its application to incompatible substructures,? Computers and Structures , 51, pp. 459?473. [47] Fonseca P. D., 2000, ?Simulation and Optimization of the Dynamic Behavior of Mechatronics Systems,?, Ph.D. Thesis, Katholieke Universiteit Leuven. [48] Bishop, R, Johnson D., 1960, The mechanics of vibration, Cambridge University Press, Cambridge. [49] Schmitz T. L., and Smith K. S., 2009, Machining Dynamics - Frequency Response to Improved Productivity, Springer. [50] Liu G. R., 2002, Mesh Free Methods: Moving Beyond the Finite Element Method, CRC Press. [51] Tlusty, J, Ziegert, J, Ridgeway S., 1999, ?Fundamental Comparisons of the Use of Serial Parallel Kinematics for Machine Tools,? CIRP Annals - Manufacturing Technology, 48(1), pp. 351?356. [52] Tosatti L. M., Bianchi G., Fassi I., Boer C. R., and Jovane F., 1998, ?An Integrated Methodology for the Design of Parallel Kinematic Machines,? CIRP Annals - Manufacturing Technology1, 47(1), pp. 341?345. [53] Oiwa T., 2010, ?Precision Mechanisms Based on Parallel Kinematics,? International Journal of Automation Technology, 4(4), pp. 326?337. [54] Zhang D., Wang L., and Lang S., 2005, ?Parallel Kinematic Machines: Design, Analysis, and Simualtion in an Integrated Virtual Environment,? ASME Journal of Mechanical Design, 127, pp. 580?588. [55] Munzinger, C, Krausse, M, Kipfmuller M., 2010, ?Simulation of parallel kinematic machine tools with minimal effort,? Production Engineering Research and Development. 125 [56] Wang, X, Mills J., 2006, ?Dynamic modeling of a flexible-link planar parallel platform using a substructuring approach,? Mechanism and Machine Theory, 41, pp. 671?687. [57] Henninger, C, and Eberhard P., 2009, ?Computation of stability diagrams for milling process with parallel kinematic machine tools,? Proc. IMechE Part I: J. Systems and Control Engineering, 223, pp. 117?129. [58] Noor A. K., 1994, ?Recent advances and applications of reduction methods,? Applied Mechanics Review, 47(5), pp. 125?146. [59] Qu Z. Q., 2004, Model Order Reduction Techniques With Applications in Finite Element Analysis, Springer-Verlag London. [60] Craig Jr, R R K. A. J., 2006, Fundamentals of Structural Dynamics, John Wiley ans Sons, INC. [61] Moore B., 1981, ?Principal component analysis in linear systems: Controllability, observability, and model reduction,? IEEE Transaction on Automatic Control, 26(1), pp. 17?32. [62] Wilcox, K, Peraire J., 2002, ?Balanced model reduction via the proper orthogonal decomposition,? AIAA Journal, 40(11), pp. 2323?2330. [63] Gawronski W., 2004, Advanced Structural Dynamics and Active Control of Structures, Springer. [64] Benner P., Bonin T., Fassbender H., Saak J., Soppa A., and Zaeh M. F., 2009, ?Novel Model Reduction Techniques for Control of Machine Tools,? ANSYS Conference & 27th CADFEM Users Meeting, Leipzig, Germany. [65] Guyan R. J., 1965, ?Reduction of stiffness and mass matrices,? AIAA Journal, 3, p. 380. [66] O?Callahan 1989 J. C., 1989, ?A procedure for an improved reduced system (IRS) model,? 7th IMAC, Las Vegas, pp. 17?21. [67] Friswell M. I., Garvey S. D., and Penny J. E. T., 1995, ?Model reduction using dynamic and iterated IRS techniques,? Journal of Sound and Vibration, 186, pp. 311?323. [68] O?Callahan, J, Avitablie, O, Riemer R., 1989, ?System Equivalent Reduction Expansion Process (SEREP),? Seventh IMAC, Las Vegas. 126 [69] Koutsovasilis P., and Beitelschmidt M., 2010, ?Model order reduction of finite element models: improved component mode synthesis,? Mathematical and Computer Modelling of Dynamical Systems, 16, pp. 57?73. [70] Givoli D., Barbone P. E., and Patlashenko I., 2004, ?Which are the important modes of a subsystem?,? International Journal for Numerical Methods in Engineering , 59, pp. 1657?1678. [71] Litz L., 1981, ?Order reduction of linear state-space models via optimal approximation of the non-dominant modes,? Large Scale Systems , 2, pp. 171?184. [72] Feliachi A., 1990, ?Identification of critical modes in power systems,? IEEE Transaction on Power Systems, 5, pp. 783?787. [73] Shabana A., and Wehage A R., 1983, ?Variable Degree-of-Freedom Component Mode Analysis of Inertia Flexible Mechanical Systems,? ASME Journal of Mechanisms, Transmissions and Automation in Design, 105, pp. 371?378. [74] Tayeb S., and Givoli D., 2011, ?Optimal modal reduction of dynamic subsystems: Extensions and improvements,? International Journal for Numerical Methods in Engineering , 85, pp. 1?30. [75] Park K. C., and Park Y. H., 2004, ?Partitioned Component Mode Synthesis via a Flexibility Approach,? AIAA Journal, 42(6), pp. 1236?1245. [76] Jakobsson H., Bengzon F., and Larson M., 2011, ?Adaptive component mode synthesis in linear elasticity,? International Journal for Numerical Methods in Engineering , 86, pp. 829?844. [77] Yuan Lin C., Pin Hung J., and Liang Lo T., 2010, ?Effect of preload of linear guides on dynamic characteristics of a vertical column?spindle system,? International Journal of Machine Tools and Manufacture, 50(8), pp. 741?746. [78] Rivin E., 1999, Stiffness and Damping in Mechanical Design, Marcel Dekker. [79] Pal, D, Basu S., 1973, ?Surface preparation of machine tool slideways and its influence on the contact stiffness,? Tribology, pp. 173?177. [80] Mi L., Yin G., Sun M., and Wang X., 2012, ?Effects of preloads on joints on dynamic stiffness of a whole machine tool structure,? Journal of Mechanical Science and Technology, 26(2), pp. 494?508. 127 [81] Dhupia J. S., Powalka B., Ulsoy a. G., and Katz R., 2007, ?Effect of a Nonlinear Joint on the Dynamic Performance of a Machine Tool,? Journal of Manufacturing Science and Engineering, 129(5), p. 943. [82] Dhupia, J S, Powalika, B, Ulsoy A G, Katz R., 2007, ?Experimental Identification of the Nonlinear Parameters of an Industrial Translational Guide for Machine Performance Evaluation,? Control, Journal of Vibration, 14, pp. 645?668. [83] Furukawa, Y, Moronuki N., 1987, ?Contact Deformation of a Machine Tool Slideway and Its Effect on Machining Accuracy,? JSME International Journal, 30(263), pp. 868?874. [84] Guo L., Zhang H., Ye P., and Zhao T., 2010, ?On Obtaining Machine Tool Joints Stiffness by Integrated Modal Analysis,? International Conference on Mechanic Automation and Control Engineering (MACE), pp. 2661?2664. [85] Kim J., Yoon J., and Kang B. S., 2007, ?Finite Element Analysis and Modeling of Structure with Bolted Joints,? Applied Mathematical Modeling, 31, pp. 895?911. [86] Natke H. G., 1988, ?Updating computational models in the frequency domain based on measured data: a survey,? Probabilistic Engineering Mechanics, 3(1), pp. 28?35. [87] Baruch M., 1978, ?Optimization Procedure to Correct Stiffness and Flexibility Matrices using Vibration Data,? AIAA Journal, 16(11), pp. 1208?1210. [88] Friswell M. I., Inman D. J., and Pilkey D. F., 1998, ?The Direct Updating of Damping and Stiffness Matrices,? AIAA Journal, 36(3), pp. 491?493. [89] Heylen W., Lammens S., and Sas P., 1997, Modal Analysis Theory and Testing, Katholieke Universiteit Leuven. [90] Crawley, E F, O?Donnell K. J., 1987, ?Force-State Mapping Identification of Nonlinear Joints,? AIAA Journal, 25(7). [91] Jalali H., Ahmadian H., and Mottershead J. E., 2007, ?Identification of nonlinear bolted lap-joint parameters by force-state mapping,? International Journal of Solids and Structures, 44(25-26), pp. 8087?8105. [92] Schmitz T. L., and Donaldson R., 2000, ?Predicting high-speed machining dynamics by substructure analysis,? CIRP Annals - Manufacturing Technology, 49, pp. 889?896. 128 [93] Park S. S., Altintas Y., and Movahhedy M., 2003, ?Receptance Coupling for end mills,? International Journal of Machine Tools and Manufacture, 43, pp. 889?896. [94] Erturk A., Ozguven H. N., and Budak E., 2006, ?Analytical modeling of spindle-tool dynamics on machine tools using Timoshenko beam model and receptance coupling for the prediction of tool point FRF,? International Journal of Machine Tools and Manufacture, 46(1901-1912). [95] Kumar U., and Schmitz T. L., 2012, ?Spindle Dynamics Identification for Receptance Coupling Substructure Analysis,? Precision Engineering, 36(3), pp. 435?443. [96] Park S., and Chae J., 2008, ?Joint Identificationof Modular Tools Using a Novel Receptance Coupling Method,? International Journal of Advanced Manufacturing Technology, 35, pp. 1251?1262. [97] Altintas Y., 2000, Manufacturing Automation: Metal Cutting Mechanics, Machine Tool Vibrations, and CNC Design, Cambridge University Press. [98] Altintas Y., and Weck M., 2004, ?Chatter Stability in Metal Cutting and Grinding,? Annals of the CIRP, 53(2), pp. 619?642. [99] Brecher C., Esser M., and Witt S., 2009, ?Interaction of Manufacturing Process and Machine Tool,? CIRP Annals, 58, pp. 588?607. [100] Park, H, Park, Y, Liang S., 2011, ?Multi-procedure design optimization and analysis of mesoscale machine tools,? International Journal of Advanced Manufacturing Technology, 56, pp. 1?12. [101] Shamoto E., and Akazawa K., 2009, ?Analytical prediction of chatter stability in ball end milling with tool inclination,? CIRP Annals - Manufacturing Technology, 58(1), pp. 351?354. [102] Shamoto E., Fujimaki S., Sencer B., Suzuki N., Kato T., and Hino R., 2012, ?A novel tool path/posture optimization concept to avoid chatter vibration in machining ? Proposed concept and its verification in turning,? CIRP Annals - Manufacturing Technology, 61(1), pp. 331?334. [103] Kilic Z. M., and Altintas Y., 2011, ?Effect of Unmatched Spindle Dynamics on the Material Removal Rates in Milling Operations,? The Proceedings of MTTRF 2011 Meeting, Chicago, Illinois, U.S.A, pp. 19?24. [104] ?ANSYS V12, 2009, Documentation for ANSYS V12.? 129 [105] Quy, N D, Matzenmiller A., 2008, ?A solid-shell element with enhanced assumed strains for higher order shear deformations in laminates,? Tech. Mech., 28, pp. 334?355. [106] Cook R. D., Malkus D. S., Plesha M. E., and Witt R. J., 2001, Concepts and Applications of Finite Element Analysis, John Wiley ans Sons, INC. [107] Cao Y., 2006, ?Modeling of High Speed Machine Tool Spindle Systems,?, Ph.D. Thesis, University of British Columbia. [108] Ewins D. J., 2000, Modal Testing: Theory, Practice, and Application, Research Studies Press. [109] Chen G., 2001, ?FE Model Validation for Structural Dynamics,? Imperial College of Science, Technology & Medicine, Ph.D. Thesis, University of London. [110] Abel, J, Shephard M., 1979, ?An algorithm for multipoint constraints in finite element analysis,? International Journal for Numerical Methods in Engineering, 14(3), pp. 464?467. [111] Heirman G., and Desmet W., 2010, ?Interface reduction of flexible bodies for efficient modeling of body flexibility in multibody dynamics,? Multibody System Dynamics , 24, pp. 219?234. [112] Liu G. R., and Quek S. S., 2003, The Finite Element Method: A Practical Course, Butterworth-Heinemann. [113] Felippa C., 1978, ?Iterative Procedures for Improving Penalty Function Soultions of Algebraic Systems,? International Journal for Numerical Methods in Engineering, 12, pp. 821?836. [114] Nour-Omid, B, Wriggers P., 1987, ?A Note on the Optimum Choice for Penalty Parameters,? Communications in Applied Numerical Methods, 3, pp. 581?585. [115] Jones A. B., 1960, ?A General Theory for Elastically Constrained Ball and Radial Roller Bearings Under Arbitary Load and Speed Conditions,? ASME Journal of Basic Engineering, pp. 309?320. [116] Jorgensen, B R and Shin Y. C., 1998, ?Dynamics of Spindle-Bearing Systems at High Speeds Including Cutting Load Effects,? ASME Journal of Manufacturing Science and Engineering, 120(2), pp. 387?394. [117] ?THK Global (www.thk.com/).? 130 [118] ?CUTPRO V9.3, Advanced Machining Simulation Software (www.malinc.com).? [119] Orfaidis S. J., 1996, Introduction to Signal Processing, Prentice-Hall, New Jersy. [120] Schmitz T. L., and Duncan G., 2005, ?Three-component receptance coupling substructure analysis for tool point dynamics prediction,? ASME Journal of Manufacturing Science and Engineering, 127, pp. 781?790. [121] Ibrahim R. A., and Pettit C. L., 2005, ?Uncertainties and dynamic problems of bolted joints and other fasteners,? Journal of Sound and Vibration, 279. [122] Altintas Y., and Budak E., 1995, ?Analytical prediction of stability lobes in milling,? Annals of CIRP, 44, pp. 357?362. [123] Fleischer J., Munzinger C., and Trondle M., 2008, ?Simulation and Optimization of Complete Mechanical Behaviour of Machine Tools,? Production Engineering Research and Development, 2, pp. 85?90. [124] ?HyperMesh V11, Altair HyperWorks.? [125] Neugebauer, R, Drossel, W G, Ihlenfeldt, S, Harzbecker C., 2009, ?Design method for machine tools with bionic inspired kinematics,? CIRP Annals - Manufacturing Technology, 58(1), pp. 371?374. [126] Schroeder, T, Krabbes, M N. R., 2007, ?Reactive trajectory splitting function for machine tools with hierarchical drive structures,? The International Journal of Advanced Manufacturing Technology, 33, pp. 988?993. [127] Hohlfeld J., Th?mmler R., Andersen O., and Jehring U., 2011, ?Hochd?mpfende, zellulare Konstruktionswerkstoffe f?r den Werkzeugmaschinenbau,? Metall, 65, pp. 32?36.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Position-dependent dynamics and stability of machine...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Position-dependent dynamics and stability of machine tools Law, Mohit 2013
pdf
Page Metadata
Item Metadata
Title | Position-dependent dynamics and stability of machine tools |
Creator |
Law, Mohit |
Publisher | University of British Columbia |
Date Issued | 2013 |
Description | Machine tool’s productivity and ability to produce a component of the required quality is directly influenced by its dynamic stiffness at the tool center point. Lack of dynamic stiffness may lead to unstable regenerative chatter vibrations which are detrimental to the performance. The chatter vibrations are influenced by the changing structural dynamics of the machine as the tool moves along the tool path, resulting in position-varying machining stability of the system. Evaluation of these varying dynamics at the design stage is a complex process, often involving the use of large order finite element (FE) models. Complexity and computational costs associated with such FE models limit the analyses to one or two design concepts and at only a few discrete positions. To facilitate rapid exploration of several design alternatives and to evaluate and optimize each of their position-dependent dynamic behavior, a generalized bottom-up reduced model substructural synthesis approach is proposed in this thesis. An improved variant of the component mode synthesis method is developed and demonstrated to represent higher order dynamics of each of the machine tool components while reducing the computational cost. Reduced substructures with position-invariant response are synthesized at their contacting interfaces using novel adaptations of constraint formulations to yield position-dependent response. The generalized formulation is used to evaluate the position-dependent behavior of two separate machine tools: one with a serial kinematic configuration, and another with hybrid serial-parallel kinematics. The reduced machine model is verified against full order models and is also validated against measurements by including joint characteristics in the model. The effects of position and feed-direction-dependent compliances on machining stability are investigated by using a novel position and feed-direction-dependent-process-stability performance criterion that evaluates the productivity of machine tools in its entire work volume. Parameters limiting the target productivity levels are identified and modified; and, the complete dynamics are rapidly re-analyzed using the developed models. Optimal design modifications are shown to increase productivity by ~35%. The proposed methods in this thesis enable efficient simulation of structural dynamics, stability assessment as well as interactions of the CNC and cutting process with the machine tool structure in a virtual environment. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2013-08-13 |
Provider | Vancouver : University of British Columbia Library |
Rights | Attribution 2.5 Canada |
DOI | 10.14288/1.0074067 |
URI | http://hdl.handle.net/2429/44804 |
Degree |
Doctor of Philosophy - PhD |
Program |
Mechanical Engineering |
Affiliation |
Applied Science, Faculty of Mechanical Engineering, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2013-11 |
Campus |
UBCV |
Scholarly Level | Graduate |
Rights URI | http://creativecommons.org/licenses/by/2.5/ca/ |
AggregatedSourceRepository | DSpace |
Download
- Media
- 24-ubc_2013_fall_law_mohit.pdf [ 5.98MB ]
- Metadata
- JSON: 24-1.0074067.json
- JSON-LD: 24-1.0074067-ld.json
- RDF/XML (Pretty): 24-1.0074067-rdf.xml
- RDF/JSON: 24-1.0074067-rdf.json
- Turtle: 24-1.0074067-turtle.txt
- N-Triples: 24-1.0074067-rdf-ntriples.txt
- Original Record: 24-1.0074067-source.json
- Full Text
- 24-1.0074067-fulltext.txt
- Citation
- 24-1.0074067.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.24.1-0074067/manifest