SMOOTH TRAJECTORY GENERATION AND PRECISION CONTROL OF 5-AXIS CNC MACHINE TOOLS by BURAK SENCER BSc, Istanbul Technical University, Istanbul, Turkey, 2003 MASc, University of British Columbia, Vancouver, Canada, 2005 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) October 2009 © Burak Sencer, 2009 ii Abstract This thesis presents experimentally verified optimal feedrate generation and high performance precision control algorithms developed for 5-axis machine tools. A feedrate scheduling algorithm has been introduced to minimize the cycle time for 5-axis machining of curved tool-paths. The variation of the feed along the tool-path is expressed in a cubic B-spline form as a function of the arc displacement. The velocity, acceleration and jerk limits of the five axis drives are considered in finding the most optimal feed along the tool-path to ensure smooth and linear operation of the servo drives with minimal tracking error. Improvement in the productivity and linear operation of the drives are demonstrated through 5-axis experiments. In an effort to design an accurate contour controller, analytical models are devel- oped to estimate the contour errors during simultaneous 5-axis machining. Two types of contouring errors are defined by considering the normal deviation of tool tip from the reference path, and the normal deviation of the tool axis orientation from the reference orientation trajectory. A novel multi-input-multi-output sliding mode controller is introduced to directly minimize the tool tip and tool orientation errors, i.e. the contour- ing errors, along the 5-axis tool-paths. The stability of the control scheme is proven analytically, and the effectiveness of this new control strategy has been demonstrated experimentally. An identification technique for identifying the closed loop transfer function of ma- chine tool feed drives has been introduced. The drive system is identified in closed loop, including the feed drive mechanism, motor amplifier, and the control law. A short Numerical Control Program is used for exciting the axis dynamics without interfering iii with the servo control loop. A generalized drive model is utilized to capture the key dynamics of the drive systems, while guaranteeing the stability of the identified model dynamics by solving a constrained optimization problem. Methods developed in this thesis have been evaluated on a table tilting 5-axis ma- chining center. Their application to other 5-axis machines would require modeling of the kinematic chain and the drive dynamics to be considered in the control law design and trajectory generation. iv Table of Contents Abstract ............................................................................................................................ ii Table of Contents ........................................................................................................... iiv  List of Tables ................................................................................................................. viii  List of Figures ................................................................................................................. ix  Nomenclature ................................................................................................................ xiii  Acknowledgments ......................................................................................................... xvi  1. Introduction ................................................................................................................. 1  2. Literature Review ........................................................................................................ 7  2.1  Introduction ........................................................................................................... 7  2.2  Feed Optimization for Multi-Axis Machine Tools ................................................ 8  2.3  Modeling of Contouring Errors in Machine Tools .............................................. 10  2.4  Contouring Controller Design ............................................................................. 14  2.5  Identification of Feed Drive Dynamics ............................................................... 17  3. Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints ....... 20  v 3.1  Introduction ......................................................................................................... 20  3.2  General Kinematic Model of 5-Axis Machine Tools .......................................... 22  3.2.1  Forward Kinematics Modeling ........................................................................ 25  3.2.2  Inverse Kinematics Solution ............................................................................ 29  3.3  Minimum Time Feed optimization Problem ....................................................... 31  3.4  Kinematic Constraints ......................................................................................... 36  3.4.1  Velocity Constraints ........................................................................................ 36  3.4.2  Acceleration Constraints .................................................................................. 37  3.4.3  Jerk Constraints ............................................................................................... 41  3.5  Solution to Minimum Time Feed Optimization Problem .................................... 44  3.6  Simulations and Experimental Results ................................................................ 47  3.7  Conclusion ........................................................................................................... 53  4. Modeling of Contouring Error for 5-Axis Machine Tools ..................................... 55  4.1  Introduction ......................................................................................................... 55  4.2  Tracking Error Dynamics on the Workpiece ....................................................... 57  4.3  Tool Tip Contour Error Estimation ..................................................................... 62  4.4  Tool Orientation Contour Error Estimation ......................................................... 65  4.5  Implementation and Experimental Results .......................................................... 68  vi 4.6  Conclusion ........................................................................................................... 75  5. Precision Contouring Controller Design for 5-Axis Machine Tools ..................... 76  5.1  Introduction ......................................................................................................... 76  5.2  Design of the Sliding Mode Controller for Tool Tip Contour Errors .................. 77  5.3  Design of Sliding Mode Controller for Tool Orientation Contour Errors ........... 84  5.4  Implementation and Experimental Results .......................................................... 92  5.5  Conclusion ......................................................................................................... 101  6. Closed Loop Identification of Feed Drive Systems ............................................... 103  6.1  Introduction ....................................................................................................... 103  6.2  Closed Loop Axis Dynamics ............................................................................. 104  6.3  Least Squares (LS) Parameter Identification ..................................................... 107  6.4  Constrained Parameter Identification Using the Particle Swarm Optimization (PSO) Method ................................................................................................................ 110  6.4.1  The Principle of PSO ..................................................................................... 112  6.4.2  PSO Solution ................................................................................................. 113  6.4.3  PSO Search Parameter Tuning ...................................................................... 116  6.5  Simulation and Experimental Results ................................................................ 118  6.5.1  Validation on a Virtual Drive ........................................................................ 119  vii 6.5.2  Validation on a Ball Screw Drive .................................................................. 124  6.5.3  Validation on a Virtual 5-Axis Machine Tool ............................................... 127  6.6  Conclusion ......................................................................................................... 127  7. Conclusions ............................................................................................................... 129  7.1  Conclusions ....................................................................................................... 129  7.2  Future Research Directions................................................................................ 131  Bibliography ................................................................................................................. 133  viii List of Tables Table 3.1: Drive limits used in the feed optimization. .................................................... 48  Table 4.1: Rigid body drive parameters. ......................................................................... 68  Table 4.2: Five-axis SMC control parameters (SISO scheme)........................................ 68  Table 5.1: Drive parameters. ........................................................................................... 93  Table 5.2: 5-Axis SMC control parameters (SISO scheme). .......................................... 93  Table 5.3: 5-axis contour controller parameters. ............................................................. 95  Table 5.4 Performance of the tool tip contour controller at different feedrates. ............. 98  Table 5.5: The effect of directional bandwidths on the tool tip controller’s performance. ......................................................................................................................................... 98  Table 6.1: Actual and estimated closed loop parameters for P-PI controlled virtual drive(simulation). .......................................................................................................... 119  Table 6.2: Actual and estimated closed loop parameters for PID controlled virtual drive(simulation). .......................................................................................................... 123  Table 6.3: Actual and estimated drive parameters of a P-PI controlled virtual 5-axis machine tool. ................................................................................................................. 126  ix List of Figures Figure 1.1: Surface finish errors during 5-axis machining. ............................................... 1  Figure 1.2: Application of high speed 5-axis machining in aerospace parts, such as impeller and compressor blades......................................................................................... 3  Figure 1.3: Application of high speed 5-axis machining in automotive and aerospace parts, dies and molds. ........................................................................................................ 6  Figure 2.1: Definition of contour errors in two-axis Cartesian machining: a) contour error (linear tool-path), b) contour error(curved tool-path). ............................................ 11  Figure 2.2: Definition of contour errors in 5-axis machining: a) tool tip contour error, b) tool orientation contour error. .......................................................................................... 13  Figure 2.3: General feed drive model. ............................................................................. 18  Figure 3.1: Tool-path and motion generation for 5-axis machine tools: a) tool motion in workpiece coordinates, b) motion generation on 5-axis machine tools. .......................... 21  Figure 3.2: Coordinate systems for modeling 5-axis machine tool kinematics. .............. 25  Figure 3.3: Tangential feed, acceleration and jerk profiles. ............................................ 32  Figure 3.4: Acceleration feasible operating region.......................................................... 39  Figure 3.5: Jerk feasible operating region. ...................................................................... 43  Figure 3.6: Feed optimization along the tool-path........................................................... 45  Figure 3.7: Windowing technique: a) forward windowing (start), b) forward windowing (shifted), c)forward windowing (finalized), c) forward winding (extended). ................ 45  Figure 3.8: Five-axis tool-path. ....................................................................................... 46  Figure 3.9: Axis reference displacement. ........................................................................ 48  x Figure 3.10: Axis limits reflected on the tangential feed: a) feedrate limited by the axis velocity constraints, b) feedrate limited by the axis acceleration constraints, c) feedrate limited by the axis jerk constraints. ................................................................................. 50  Figure 3.11: Optimized and un-optimized feed profiles. ................................................. 51  Figure 3.12: Optimized and un-optimized axis profiles. ................................................. 52  Figure 3.13: Tracking errors of the axes. ......................................................................... 53  Figure 4.1: Definition of contour errors in 5-axis machining: a) tool motion in workpiece coordinates (P-system), b) contour error definition (orientation and tool tip contour errors). .............................................................................................................................. 56  Figure 4.2: Linear rigid body feed drive dynamics. ........................................................ 60  Figure 4.3: Modeling of tool tip contour errors in 5-axis machining. ............................. 62  Figure 4.4: Modeling of tool orientation contour errors in 5-axis machining. ................ 66  Figure 4.5: Test Tool-path: a) reference tool tip ( ), , ,, , Tref x ref y ref z refP P P P⎡ ⎤= ⎣ ⎦ position, b) reference tool axis orientation ( ), , ,, , Tref i ref j ref k refO O O O⎡ ⎤= ⎣ ⎦ .......................................... 69  Figure 4.6: Comparison of tool tip contouring error estimations (experimental): a) comparison of true and estimated tool tip contour errors, b) discrepancy of proposed estimation from true contour errors, c) delay time .......................................................... 71  Figure 4.7: Comparison of tool orientation contour error estimations: a) A axis errors, b) C axis errors, c) discrepancy in the proposed contour error estimation in A axis, d) discrepancy in the proposed contour error estimation in C axis. ..................................... 73  Figure 5.1: Tool tip contour error estimation. ................................................................. 78  Figure 5.2: Tool tip contour controller block diagram. ................................................... 83  xi Figure 5.3: Tool orientation contour controller block diagram. ...................................... 87  Figure 5.4: Complete 5-axis contour controller block diagram. ...................................... 92  Figure 5.5: Five-axis tool-path. ....................................................................................... 93  Figure 5.6: Reference axis trajectories: reference Cartesian axis positions, b) reference rotary axis positions, c) reference Cartesian axis velocities, d) reference rotary axis velocities. ......................................................................................................................... 94  Figure 5.7: Contouring performance of the SISO controller (experimental): a) tool tip contour errors, b) X axis tracking, c) Y axis tracking, d) Z axis tracking, e) A axis tracking, f) C axis tracking. ............................................................................................. 95  Figure 5.8: Proposed tool tip contour controller performance (experimental) ): a) tool tip contour errors, b) X axis tracking, c) Y axis tracking, d) Z axis tracking, e) A axis tracking, f) C axis tracking. ............................................................................................. 97  Figure 5.9: Effect of controller parameters.................................................................... 100  Figure 5.10: Orientation contour controller performance at F=50[mm/sec] (experimental): a) orientation contour error, b) mean orientation contour error trend. . 101  Figure 6.1: Rigid body axis dynamics and common CNC drive control structures: a) rigid body drive model, b) P-PI control, c) PID control. ............................................... 104  Figure 6.2: General closed loop dynamics model. ........................................................ 107  Figure 6.3: Identification scheme. ................................................................................. 108  Figure 6.4: Allowable locations for; a) real, b) complex poles. .................................... 112  Figure 6.5: Search region for the Particle Swarm Optimization algorithm. .................. 113  Figure 6.6: Reference identification trajectory send to the drives. ................................ 118  xii Figure 6.7: True and estimated FRFs for P-PI controlled drive system (simulation): a) tracking transfer function, b) disturbance transfer function. ......................................... 120  Figure 6.8: Parameter convergence of the PSO algorithm for P-PI controlled drive system. ........................................................................................................................... 121  Figure 6.9: Predicted and actual tracking performance of the P-PI controlled virtual feed drive system (simulation). ............................................................................................. 122  Figure 6.10: True and estimated FRFs for PID controlled drive system (simulation): a) tracking transfer function, b) disturbance transfer function. ......................................... 123  Figure 6.11: Predicted and actual tracking performance of the PID controlled virtual drive system (simulation). ............................................................................................. 124  Figure 6.12: Predicted and experimentally verified tracking performance of P-PI (a) and SMC (b) controlled ball screw drive: a) P-PI control results, b) SMC control results. . 125  Figure 6.13: Predicted and actual tracking Performance of the P-PI controlled virtual 5- axis machine tool (Simulation). ..................................................................................... 127  xiii Nomenclature , , , , , ( )X Y Z A C B Five-axis CNC machine tool axis variables , , , , , ( )x y z a c b Five-axis CNC machine tool axis position, x,y,z [mm], a,c,(b) [rad]. s Path displacement along the spline tool-path [mm] , ,s s s   Tangential velocity (feed), acceleration and jerk profiles as a function of the path displacement along the path. if Feed at B-spline control points [ ]Tx y zP P P P= Tool tip position vector and scalar components [mm] [ ]Ti j kO O O O= Unit tool axis orientation vector and its principal components in x,y,z directions [ , , , , , ]Tx y z i j kR P P P O O O= Tool motion interpolation vector function and its scalar components (in the workpiece coordinates) [ , , , , , ( )]Tx y z a c b=Q Axis position interpolation vector function and its five scalar components 5 1×Q , 5 1×Q , 5 1×Q Vectors containing the axis velocity, acceleration and jerk profiles [ , , , , ] , [ , , , , ] , [ , , , , ] T s s s s s T ss ss ss ss ss T sss sss sss sss sss x y z a c x y z a c x y z a c = = = s ss sss Q Q Q Vectors containing the first, second and third derivatives of axis positions with respect to s xiv u Spline parameter u Knot vector of the B-spline tool-path. [ , , , , ]Tx y zP P P a c=x Tool pose vector , , ,[ , , , , ] T ref x ref y ref z ref ref refP P P a c=x Reference tool pose vector [ , , , , ]TPx Py Pz a ce e e e e=e Tool pose tracking error vector in workpiece frame [ , , , , ]Tt n b a ce e e e e=Fe Tracking error vector in the Frenet Frame Re Rotary drive tracking error vector oe Tool orientation tracking vector in spherical coordinates F Frenet frame transformation matrix J , oJ Jacobian matrix, orientation jacobian ε Contour error vector in cartesian coordinates oε Contour error vector in spherical coordinates Rε Contour error vector in rotary drive coordinates ω Angular velocity of the tool axis v Velocity vector of the drives I P DC ,C ,C Integral, proportional and derivative gains S Sliding surface sK Diagonal feedback gain matrix V Lyapunov function xv d̂ Estimated disturbance κ Diagonal gain matrix , nζ ω Damping ratio and natural frequency W Diagonal weight matrix ,1 ,2 ,[ , , ] T i i i i Dσ σ σ=σ … Optimization particles ,1 ,2 ,[ , , ] T i i i i Dv v v=v … Velocity vector of the particles CNC Computer numerical control CAD Computer aided drawing CAM Computer aided manufacturing B-spline Bezier spline NURBS Non-uniform rational B-spline PI Proportional integral controller PID Proportional integral derivative controller SMC Sliding mode controller SISO Single-input-single-output MIMO Multi-input-multi-output PSO Particle swarm optimization LS Least squares xvi Acknowledgments I would like to express my sincere gratitude to my research supervisor Dr. Yusuf Altintas for his valuable instruction, guidance, and support which he has provided throughout my graduate studies at the Manufacturing Automation Laboratory. I wish to thank my colleagues in the Manufacturing Automation Laboratory for sharing their knowledge and experience with me. I particularly feel luckily to gain such an exceptional friend, Doruk. He was with me especially during the tough times and made my life colourful and enjoyable. Finally, I am deeply grateful to my beloved family, my mother Serap and my father Necip for their love, encouragement, endless support and patience they have provided me throughout my entire life and during my abroad studies. I dedicate this thesis to my dear mother, Serap. 1 1. Introduction The ongoing and increasing demands for productivity and quality improvements in manufacturing processes have focused the attention of the machine tool community on high speed machining (HSM). HSM requires multi-axis machine tools to have high feedrate and acceleration capabilities and still maintain micron positioning accuracy. At the same time, part geometries have become increasingly complex. Free-form sculptured surfaces are widely used in die and mold machining, especially in automotive and aerospace industries where manufacturers are requiring high speed machining of parts. High speed machining of free-form surfaces with continuously changing curvatures and deep cavities can only be performed on the recently developed 5-axis computer numeri- cal controlled (CNC) machines. The multi-axis machine tools have a non-orthogonal kinematic structure, and accommodate the continuous change of the tool orientation that is necessary to follow complex surface contours. Achieving the desired high speed feed motion accurately on 5-axis machine tools requires a significant amount of effort to be Figure 1.1: Surface finish errors during 5-axis machining. Chapter 1: Introduction 2 dedicated to motion planning and control law development. A surface finish error during a typical 5-axis machining operation is illustrated in Figure 1.1. The goal is to keep the actual tool position on the desired tool-path at all times, which satisfies the tolerance of the part being machined. In practice, the tool motion deviates from the desired trajectory due to the limited bandwidth of the servo drives, as well as disturbances originating from cutting forces and guideway friction. Each axis tracks its reference position command individually with a tracking er- ror ( , ,....)x ye e . As a result, the actual toolpath exhibits a geometric deviation from the reference path, represent the tool tip error, ε , and tool orientation error, θε , in Figure 1.1. The most commonly used method used to decrease the contour error is to increase the axis servo bandwidth. However, the control bandwidth can only be increased until the unmodeled dynamics in the drive system start to cause instabilities. When tracking a line segment using the Cartesian axis at constant speed, the contour error remains constant and linearly proportional to the axis tracking errors. Therefore, any improve- ment in the tracking bandwidth pushes the tool closer to the reference tool-path. On the other hand, in case of tracking a non-linear spline trajectory on a 5-axis machine tool, the linear relationship between the contour and the individual joint tracking errors is nonli- near. Because of the non-orthogonal kinematic configuration, the contribution of each axis’ tracking error varies along the tool-path, and small tracking errors from rotary drives may induce large tool deviations. Hence, in order to design an accurate contour controller for 5-axis machine tools, the contribution of individual axis tracking errors to the overall contour error must be identified, and the incorporation of the 5-axis machine tool kinematics into the servo controller design becomes essential. Chapter 1: Introduction 3 Spline tool-paths have recently gained importance for delivering the desired feed motion smoothly and accurately during complex sculptured surface machining processes (Figure 1.2). The feed profile along a spline tool-path should satisfy particular require- ments for HSM including: maximum smoothness for better tracking performance, and minimum cycle time to increase productivity. For instance, selecting the highest cruise velocity along the toolpath is ideal for productivity; however, any discontinuities, or sharp changes in the spline curvature will appear in the axis kinematic profiles as high frequency spikes. As a result, tracking performance will be degraded, and, if the axis kinematic limits are violated, severe contour errors will occur. In general practice, tangential feedrate is conservatively selected based on the kinematic limits of the drives. However, this approach compromises the minimum cycle time objective. The tangential Figure 1.2: Application of high speed 5-axis machining in aerospace parts, such as impeller and compressor blades. Chapter 1: Introduction 4 feedrate should be scheduled smoothly along the tool-path. By increasing the feed at straight sections, the overall cycle time is minimized while still respecting the machine tool’s velocity, acceleration, and jerk limitations, keeping the tracking accuracy accepta- ble. Hence, scheduling the optimal feedrate profile along spline tool-paths for minimum time is not a trivial task on 5-axis machine tools. However, if realized in a practical manner, full physical advantage of the machine is employed maximizing the potential for increasing the process productivity. 5-Axis machine tools are mainly used in machining dies and molds that have smooth surfaces. During following such a smooth contour, all the axes of the machine tool are in motion. The major problems in achieving accurate high speed contour ma- chining on 5-axis machine tools, namely: smooth trajectory generation, feedrate optimization, and precision contour controller design have been tackled in this thesis. An efficient feedrate scheduling method has been developed for 5-axis machine tools executing spline tool-paths. The developed algorithms minimize the machining cycle time while adhering to the velocity, acceleration and jerk limits of the drives, which helps to preserve high speed contouring accuracy while executing an optimal feed motion. An alternative control philosophy is then introduced for improving the contour- ing accuracy of 5-axis machine tools at high speeds. Rather than tracking individual drives, the minimization of contouring error is considered as the prime objective in the control law. Two types of contour errors, namely the tool tip error and the tool orienta- tion error are defined, and analytical models are developed to estimate these errors in real time. Utilizing the machine kinematics and the contour error models, a novel multi- input-multi-output (MIMO) continuous time, integral Sliding Mode Contour controller is Chapter 1: Introduction 5 introduced. This controller is more robust against disturbances and modeling errors as compared to the widely used feedback controllers. Effectiveness of the proposed control strategy is demonstrated experimentally on the open architecture 5-axis CNC machine tool shown in Figure 1.3. Fadal VMS 2023 3-axis vertical Cartesian machine tool has been retrofitted into a 5-axis machine tool by mount- ing Nikken tilting-rotary drive system, i.e. a trunnion table. The machine table translates in Cartesian X and Y directions and carries the rotary-tilting drive assembly, as shown in Figure 1.3. The tilting axis is called A axis and introduces rotation around the X direc- tion. The spindle axis translates only in the Z direction. The rotary axis of the trunnion table is called C axis, and it introduces rotation around the Z axis of the machine tool. PWM current amplifiers are used to generate currents to the motors that deliver corres- ponding torque to 5 of the machine axes. On this experimental setup, axis amplifiers are set to operate in torque mode. Hence, an external voltage send to amplifiers generates corresponding torque with respect to the amplifier’s current gain and motor’s torque constant. Real time control of axis position is then performed as follows: Labview® real time software is used to target a Personal Computer(PC) equipped with National Instru- ments data acquisition cards to encode digital rotary encoder signal from machine tool’s axes and to generate analog control signal delivered to the amplifiers. Using Labview® real time software, mathematical computations can be performed, and various control algorithms can be implemented to run at desired sampling rates. As shown in Figure 1.3, proportional-integral-derivative (PID) or sliding mode controllers (SMC) can be pro- grammed to control position of the 5 axes in real time. Further details on drive specifications are presented in following chapters of the thesis. Chapter 1: Introduction 6 The thesis is organized as follows: Related literature is reviewed in Chapter 2, fol- lowed by general kinematic model of 5-axis machines and the feedrate scheduling method in Chapter 3. The contour error models are presented in Chapter 4 and the precision contour controller design in Chapter 5. In Chapter 6, a closed loop drive identification method is introduced, and conclusions with possible future research directions are highlighted in Chapter 7. -Workstation PC Read NC-Program -Real Time Control PC (1kHz) Decode NC-Program Generate Trajectory linear, circular, spline Position Control PID, ZPECT, SMC, Fric.Comp Ethernet - FPGA Board (40kHz) Read ADC Quadrature Encoder Decoding Write DAC PCI Bus La bv ie w R ea l T im e Encoder Feedback Motor Current Command C A Y Z Figure 1.3: Application of high speed 5-axis machining in automotive and aerospace parts, dies and molds. 7 2. Literature Review 2.1 Introduction Five-axis computer numerical controlled (CNC) machine tools have been used widely in machining parts with complex sculptured surfaces such as dies, molds and aerospace parts. The ultimate goal is to maintain the desired dimensional accuracy in high speed 5-axis machining applications. During 5-axis machining, the dimensional accuracy of the part is affected mainly by the tool positioning errors corresponding to the feed drive system tracking accuracy. The servo system must possess high closed loop bandwidth in order to accurately track the rapidly varying reference tool position com- mands. Furthermore, due to the coupled nonlinear kinematics of the 5-axis machine tool, small positioning errors of the drives may cause excessive deviation of the tool tip from the desired trajectory. Hence, design of a centralized controller can establish coordina- tion between the axes in order to achieve the desired contouring accuracy. As well, the reference tool trajectory must be smooth in order to avoid exciting the structural dynam- ics of the drive system, and thus avoid excessive contouring errors caused by motion command demands beyond the bandwidth of the axes servos. Therefore, smooth time optimal trajectory generation and accurate motion control are the two essential compo- nents for achieving contouring accuracy during 5-axis machining that are studied in this thesis. In this chapter, previous research in literature is reviewed in relation to the intended contributions in this thesis in the following areas: Feedrate optimization for multi-axis machine tools in Chapter 3, modeling of contour errors in Chapter 4, design of contour controllers in Chapter 5 and Modeling and identification of feed drives in Chapter 6. Chapter 2: Literature Review 8 2.2 Feed Optimization for Multi-Axis Machine Tools Contouring tool-paths for 5-axis machining are generated by computer aided de- sign/manufacturing (CAD/CAM) systems in cubic or non-uniform rational Bezier (NURB) spline format, which can be processed by modern CNC systems. Although a constant feed is programmed in the NC block, the velocity, acceleration and jerk profiles on all five drives may vary continuously as a function of path curvature and the inverse kinematics of the machine tool. If the velocity, acceleration and jerk of each drive exceed its physical saturation limits, the linear operation of the servo control is violated, leading to severe marks on the surface or instability of the controller. Conservative, constant feed speeds are computed by the process planner so that the axis limits are never violated; this may lead to lengthy machining times. Thus, the challenge is to generate smooth trajectories that traverse the spline tool-path in minimum time without violating the physical limitations of the machine’s drives. The identification of optimal feed profiles without exceeding the saturation limits of the actuators is a nontrivial optimization problem, and studied mainly by robotics researchers. Bobrow et al. [1] and Shin and McKay [2] were the first to solve the mini- mum time control problem. They considered the trajectory generation as a dynamic system with two states, i.e. path displacement, ( )s t , and velocity, ( )s t , and proved that the actuator torque constraints limit the velocity along the tool-path. They switched the path acceleration, ( )s t , between its maximum and minimum limits at computed switch- ing points to generate a bang-bang minimum time trajectory. Shiller and Lu [3],[4] developed a reliable search procedure to find the switching points. Other methods, such Chapter 2: Literature Review 9 as Dynamic Programming [5] or Pontryagin’s minimum principle [6] have also been applied to solve the optimal control problem more efficiently. The machine tool literature has implemented some of those approaches in the Carte- sian machining of spline tool-paths. Smith et al. [7] used the acceleration limits to generate the velocity limiting curve, and implemented a similar approach to Shiller and Lu to find the optimal feed-rate profile along spline tool-paths for a two-axis Cartesian case. Renton and Elbestawi [8] outlined a computationally efficient two-pass algorithm to solve the minimum time feed optimization problem. The algorithm works by scanning through the trajectory in the forward and reverse directions. Dong and Stori [9] have also included a series of machine tool capability constraints, and presented a two-pass structure. Unfortunately, such trajectories possess discontinuous acceleration and torque profiles, defeating the essential purpose of using a smooth tool-path and leading to inaccurate contouring in high speed machining. To achieve smooth trajectories that enhance the tracking performance, the jerk lim- its of the drives have been considered as constraints at the expense of increased computational complexity. These approaches require iterative solutions as opposed to optimal control formulations [10]. For example, Piazzi and Visioli [11] used axis jerk; Constantinescu and Croft [12] used the torque rate constraints and iteratively optimized the spline trajectories for minimum time motion of robotic manipulators. Jerk limited trajectories have also been adapted for machine tools, and computationally efficient methods have always been favored for implementation on real time CNCs. A practical feed modulation solution was proposed by Weck et al. [13] for multi-axis machine tools following quintic spline tool-paths. They considered velocity and acceleration limits of Chapter 2: Literature Review 10 the machine tool axis to determine a conservative feed for each spline segment of the tool-path. The discrete feeds are connected with jerk limited cubic acceleration trajectory profile without considering the jerk limits of the drives. Altintas and Erkorkmaz [14] later defined the displacement profile as a minimum jerk quintic spline of time along the tool-path, and optimized the travel time of each segment iteratively. They considered velocity, acceleration, and the jerk limits of the drives for a smoother motion. Although their approach is extensible to any type of multi-axis machine tool, the number of constraints to be considered increases dramatically with the number of axes. For in- stance, 14 constraints are used in the two-axis case: this increased to 32 for 5-axis machines. In this thesis, a jerk continuous feedrate optimization technique is introduced for minimizing the cycle time in 5-axis machining, while preserving the motion smoothness and tracking accuracy for traveling along spline tool-paths. The feedrate optimization technique developed for 5-axis machining has been published in [15]. 2.3 Modeling of Contouring Errors in Machine Tools Each drive on the machine tool is subject to unavoidable tracking errors between the commanded and actual positions due to the servo controller dynamics. These tracking errors are kinematically transformed to the tool tip, resulting in contouring errors be- tween the commanded tool-path and the actual path delivered by the CNC unit. In praxis, dedicated control laws of each drive try to minimize the individual axis tracking errors ( ), , ...x y ze e e . As a consequence of the nonlinear kinematics of the 5-axis machine tool, small tracking errors on rotary axes may result in larger deviations of the tool tip Chapter 2: Literature Review 11 from the reference tool-path relative to the workpiece. In addition, when one axis is subject to a disturbance, other axes will perform as if the disturbed axis is functioning normally. This results in a loss of coordination and indegradation of the overall contour- ing accuracy. The contour error ( )ε , which is defined as the orthogonal deviation from the de- sired tool-path [16] , is the actual measure of part accuracy. Figure 2.1 illustrates the relationship between contour error and axis tracking errors for simple linear and curved tool-paths. When following a linear tool-path (Figure 2.1a), the contour error ( )ε can be expressed as a simple analytical function of the axis tracking errors, ,x ye e , as: ε sin(θ) . cos(θ) x x y y x y C e C e C C = − + ⎫⎪= ⎬⎪= ⎭ (2.1) On the other hand, while machining curved tool-paths such as B-splines or NURBS, accurate computation of the contour errors becomes more challenging (Figure 2.1b). Koren et al. [17] used Equation (2.1) in estimating the contour errors in two-axis machining along simple linear tool-paths. They proposed the first cross coupling control y a) b) x ε ε1 ε2 θ Reference Position Reference Position Actual Tool Position y x De sir ed Co nto ur ex ex1 ex2 ey2 ey1 ey Figure 2.1: Definition of contour errors in two-axis Cartesian machining: a) contour error (linear tool-path), b) contour error(curved tool-path). Chapter 2: Literature Review 12 (CCC) algorithm to reduce the normal deviation of the tool tip from the path. Later, Koren and Lo [18] extended contour error estimation to design a CCC algorithm for circular tool-paths. Erkorkmaz and Altintas [19] estimated contour error along arbitrary tool-paths by computing the Euclidian distance from the desired tool-path to the actual tool position in the vicinity of a number of interpolated points, and identified the refer- ence point with the minimum deviation from the actual tool tip. The contour error vector was then defined to be normal to the vector between the identified reference point and the consecutive reference path points. As a result, the shortest vector normal to the desired contour was obtained, but the computationally intensive nature of the iterative solution proposed was undesirable for real-time control. Chiu and Tomizuka [20] developed an analytical contour error estimation procedure for two-axis machining of C2 continuous curved tool-paths. They reflected the tracking errors of the Cartesian axes on the tangent and normal directions of the tool-path, and considered the error component in the normal direction as the approximation of actual contour error. Later, Yeh and Hsu [21] defined the gains ( ) ( ),x yC t C t in Equation (2.1) to vary as a function of time. Thus, they were also able to estimate contour errors for splines and other curved paths. The contour error approximation becomes inaccurate if the error in tangent direction is high. This occurs particularly during high speed machining of paths with high curvature. Peng and Chen [22] attacked this problem and modified previous approximations with addi- tion of an extra term so that the over-cutting of the path is avoided in two-axis Cartesian motions. Apart from Cartesian contouring systems, Katz et al. [23] addressed the con- touring problem on a non-orthogonal reconfigurable three-axis machine tool. They modeled the kinematics of the machine tool and designed time varying gains for Equa- Chapter 2: Literature Review 13 tion (2.1) in order to estimate contour errors. Lin and Hu [24] used the kinematics of a 5- axis machine tool and estimated contour errors of the tool tip as the difference between the reference and actual circle radii during simple circular interpolation. Later, Lo [25] designed coupled proportional-integral-derivative (PID) controllers and showed that minimizing the tracking errors of the tool tip in the workpiece coordinates results in better contouring performance. Recent literature has mainly addressed the contouring control problem for 3-axis Cartesian machining applications, where the contour error is defined as the norm of the translational deviation of the tool tip from the reference tool-path. However in 5-axis machining, the tool translates and rotates at the same time to follow the reference orientation. Thus, two types of contour errors are proposed in this thesis. The first type considers the normal deviation of tool tip from the reference path in Cartesian coordi- nates. This is defined as the tool tip contour error. The second contouring error is the deviation of the tool axis orientation from the reference orientation trajectory in spheri- cal coordinates, called the tool orientation contour error (Figure 2.2). The contributions of the contour error models developed for the 5-axis machine tools are published in [25] and further discussed in Chapter 4. a) b) ε Actual Tool Tip PointSurface φ Actual Tool Orientation Figure 2.2: Definition of contour errors in 5-axis machining: a) tool tip contour error, b) tool orientation contour error. Chapter 2: Literature Review 14 2.4 Contouring Controller Design Two major approaches have been adopted to reduce contouring errors in machine tool literature. In the first approach, contour errors are reduced indirectly by attempting to reduce axis tracking errors. Traditional algorithms such as Proportional (P), Propor- tional-Integral (PI) and Proportional-Integral-Derivative (PID) are based on the feedback principle [27], [28]. Tomizuka [29] has developed a zero phase error tracking controller (ZPETC) that cancels the stable dynamics of the servo drive in a feed-forward fashion. The bandwidth of the overall system and thus, the tracking accuracy of the drives increases with ZPETC, provided that the drive model is accurate and does not vary with time [30], [31]. Recent efforts have been directed towards improving the bandwidth of the axes drives using sliding mode controllers that are more robust to changes in the drive dynamics. The general sliding mode controller (SMC) was first introduced by Utkin [32], and required switching around the sliding surface resulting in a disconti- nuous control law. To overcome this problem, Slotine and Li [33] proposed an adaptive sliding mode controller, which estimates and cancels various uncertainties that do not vanish at the equilibrium point. In parallel, Stepanenko et al. [34] noted that the transient response of the discontinuous SMCs could be improved by including the integral, position and derivative of tracking errors in the sliding surface design. Later, Altintas et al. [35] proposed a continuous sliding mode controller for accurate axis tracking control by considering only the rigid body dynamics of the ball screw drive. Their SMC design was based on penalizing the position and velocity errors of the axes, and demonstrated similar performance to the ZPETC, but with improved robustness to uncertainties in the drive dynamics. Chapter 2: Literature Review 15 An improved method of controlling the machine tool feed drive system is to intro- duce coupling actions in the controller to maintain coordination along the desired contour. By introducing coupling effects among multiple axes, coordinated motion is achieved by either the “equal status” or the “master-slave” approach [36]. When the dynamics are significantly different among multiple axes, the controller designed by the “equal status” approach may saturate the axes of the slower actuator. In order to over- come this drawback, the “master-slave” approach is favored. This approach assigns the slowest axis as the master of the faster drive. Su et al. [37], [38] have developed an adaptive coordination controller for position synchronization of multiple axes. They defined the synchronization error as the differential position error among multiple motion axes, and penalized this error in the feedback loop. The second major approach in reducing contour errors attempts to estimate the con- tour errors in real time and generate control action against them. Koren [17] estimated the contour error in two-axis machines as a function of axis tracking errors and linear feed direction between consecutive path points. He reduced contour errors by injecting corresponding coupling command to the servos in order to push the actual tool position towards the desired tool-path; creating the Cross Coupled Controller (CCC) mentioned above. Later, Koren and Lo [18] used time varying coupling gains to implement CCC along circular paths. When the axis controllers minimize the tracking errors along a curved path, the contour error is increased, forcing the contour controller to act against this error. As a result, it is difficult to distinguish which control element dominates the final contouring result. The analysis hindered implementation of the CCC scheme in non-orthogonal machine tools. Erkorkmaz and Altintas [19] developed a numerical Chapter 2: Literature Review 16 method to estimate the contour error for arbitrarily shaped tool-paths. They have imple- mented it in a CCC scheme together with feed-forward axis dynamics compensation and demonstrated improved contouring performance for Cartesian machining. Chiu and Tomizuka [40] used the coordinate transformation approach to directly design for the desired contour error dynamics. Using linear time varying PD regulators, the decoupled error dynamics in tangential and normal directions were stabilized. This approach demonstrated that contouring performance is improved by increasing closed loop bandwidth in the normal direction. However, the PD controller is not sufficiently robust enough and the contour error approximation may become inaccurate. Therefore, con- touring accuracy degrades in high speed machining of paths with sharp curvatures. Peng and Chen [41] defined a geometric contouring index (CI) that improves error estimation on circular tool-paths, and designed a back stepping sliding mode controller in the normal direction to introduce robustness against friction. Chen et al. [42] later designed robust controllers in polar coordinates to establish contouring control for very simple non-circular tool-paths. By utilizing the kinematics of the machine tool and contour error models presented in this thesis, a multi-input-multi-output continuous time integral sliding mode contour controller is introduced for 5-axis machine tools; this controller is robust against distur- bances and modeling errors. The effectiveness of the control strategy is compared against the conventionally used, axis based SMC and the contributions are presented in [43]. Chapter 2: Literature Review 17 2.5 Identification of Feed Drive Dynamics There has been extensive research in literature dedicated to building models of feed drive systems. One of the basic models, proposed by Koren [44], considers only linear rigid body dynamics as illustrated in Figure 2.3. In this model, the controller generates the current command, u , which is amplified by the current amplifier constant, aK , to generate the armature current, i. The armature is then converted into actuation torque, mT , via the motor constant, tK . Part of the actuation torque is consumed by disturbances ( )dT originating from nonlinear guideway friction and, possibly, the cutting process. The remaining torque actuates the drive system, and the linear rigid body model consists of an equivalent inertia, J, and the viscous friction, B, reflected by the actuator. Models have also been proposed to capture and integrate the nonlinear effect of fric- tion in the feed drive. Armstrong et al. [45] have presented a very broad survey of the available techniques used to model and compensate friction in servo systems. They have proposed a seven parameter model which adequately captures the effect of pre-sliding displacement, static friction, and the full Stribeck curve describing how friction force changes with the relative sliding velocity between two contact surfaces. Later, Lee and Tomizuka [46] have applied the effect of pre-sliding displacement, static friction, and the full Stribeck curve describing how the friction force changes with the relative sliding velocity between two surfaces in contact. Later, Erkorkmaz and Altintas [47] applied simpler versions of this model to identify the friction dynamics in CNC drives. Chapter 2: Literature Review 18 The main disadvantage of rigid body type models is that they fail to capture the ef- fect of structural vibrations. Varanasi and Nayfeh [48] and Erkorkmaz and Kamalzadeh [49] have built finite element models and constructed frequency response experiments to identify the torsional and axial vibrations of ball screw drives. Backlash and motion loss in the nut interface have been identified by Kao et al. [50] and Cuttino et al. [51], with models of varying complexity. The volumetric and thermal errors, which influence the final tool positioning accuracy, were modeled by researchers at NIST in [52]. Upon development of models to represent the key dynamic behavior of the feed drive dynamics, open and closed loop identification methods have been proposed to fit the model parameters. Altintas and Erkorkmaz [47] have developed an open loop identification scheme for feed drive systems. They send a rich identification signal to the amplifier and measured the actual position of the axis in open loop. Later, linear least squares method [53] was used to fit the inertia, and friction parameters accurately to the rigid body model. Apart from the open loop identification, closed loop identification has been favored in the industry, since it does not interrupt the machine or the process. The use of Laguerre series expansion has been widely used in order to obtain reliable data Figure 2.3: General feed drive model. Chapter 2: Literature Review 19 from noisy measurements [54], [55], [56]. Ramirez et al. [57], [58] have also developed an algebraic framework for continuous-time linear systems identification in the presence of noise. However, those identification techniques require particular identification signals to be send to the system such as step, ramp or pseudo random binary sequence (PRBS) inputs. On a real CNC, only predefined smooth trajectory commands can be used. Erkorkmaz and Wilson [59],[60],[61] developed a rapid identification method to estimate the closed loop transfer function of CNC feed drives. They used linear motion commands to generate a reference identification signal and estimated close loop drive parameters. In this thesis, the framework developed by Erkorkmaz [60] has been utilized to es- timate the closed loop dynamics of CNC feed drives. In order to guarantee the stability of the estimated model, closed loop pole locations are constrained, and a particle swarm optimization method has been developed to solve the constrained transfer function identification problem [62]. 3. Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 3.1 Introduction High speed 5-axis tool-paths for machining of dies and molds with sculptured sur- faces are generated by CAD/CAM systems as cubic or NURB splines, which can be processed by modern CNC systems. A constant cutting feed is programmed for segments of the tool-path by the process planner. However, the velocity, acceleration and jerk on all five drives may continuously vary as a function of path curvature and the inverse kinematics of the machine tool. If the velocity, acceleration or jerk of each drive exceeds its physical saturation limits, the linear operation of the servo control is violated, leading to severe marks on the surface or instability of the controller. On the other hand, choos- ing conservative speeds increases overall cycle time and cuts back the productivity severely. In this chapter, a feedrate optimization technique is introduced for minimizing the cycle time required to traverse 5-axis tool-paths without violating the physical limits of the drives. The feed profile is defined in C2 continuous cubic B-spline format as a function of the arc displacement to provide a jerk continuous motion along the tool-path. Minimum cycle time is obtained by iteratively modulating the feed control points of the B-spline feed profile to maximize path velocity while respecting process and machine tool constraints. Long toolpaths are handled efficiently by applying a moving window technique. Hence, the proposed algorithm can be used either off-line, or can be inte- grated into the CNC executive for high speed contour machining operations where limited feedrate fluctuations are not detrimental to the manufacturing process such as rough and semi-finishing of dies and molds. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 21 Five-axis machine tools exhibit non-Cartesian, nonlinear kinematics, which require mapping of trajectories planned on the workpiece into reference axis position com- mands. The same mapping must be used in transforming the physical limits of the drive system as constraints on the tangential feedrate. Therefore, a general kinematic modeling technique for widely used 5-axis machine tools is first introduced in Section 3.2. The formulation of the optimization problem is presented in Section 3.3. Kinematic con- strains are formulated in Section 3.4, and the solution is presented in Section 3.5. The algorithm is demonstrated with simulations and experiments on 5-axis flank milling of an impeller blade in Section 3.6, and the contributions [15] are subsequently summa- rized. Table Workpiece To ol R(s) = [P(s),O(s)]T (Cutter Motion) Inverse Kinematics Transformation q(s) = [x(s),y(s),z(s),A(s),C(s)]T (Machine Axis Motion) x y + + + + + z C B Xw Yw Zw s P2 P3 S3 s2 s1 S1 P1 PNp PNp-1 P=P(s) O=O(s) o2 o3 o1 ds dt PNs To ol ds dt = s : Tangential Feedrate Workpiece Coordinate System a) b) Machine Coordinate System Figure 3.1: Tool-path and motion generation for 5-axis machine tools: a) tool motion in workpiece coordinates, b) motion generation on 5-axis machine tools. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 22 3.2 General Kinematic Model of 5-Axis Machine Tools A typical 5-axis motion command, shown in Figure 3.1a, is given in the workpiece coordinate system (WCS) as a sequence of discrete positions or in spline format of continuous cutter motion along the tool-path. The tool motion is defined by the three Cartesian coordinates of its tool center point (TCP), T x y zP P P⎡ ⎤= ⎣ ⎦P , and the orienta- tion vector of the cutter axis, T i j kO O O⎡ ⎤= ⎣ ⎦O , 1=O . These coordinates are concatenated into the 6-tuple, [ ] TT x y z i j kP P P O O O⎡ ⎤= → ⎣ ⎦R P O , (3.1) in the cutter location (CL) file. Five-axis machines are composed from three Cartesian ( , , )X Y Z and two rotary ( ), or AB BC AC drives. For a given tool position and orienta- tion ( )R , the corresponding physical axis reference position command, [ ], , , , Tx y z a c=Q , (3.2) is computed from the kinematic model of the machine tool. The transformation is known as the inverse kinematics [63], which maps the spatial cutter motion from WCS into the machine tool axis (joint) system (Figure 3.1b), namely, →R Q , (3.3) and the mapping of the machine’s joint motion to the tool motion, →Q R (3.4) Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 23 is referred to as the forward kinematics transformation. Typical configurations of 5-axis machine tools can be divided into three basic types [64], [65]: (1) table-tilting (Figure 3.1b), with two rotations on the table; (2) spindle-tilting, with two rotations on the spindle, (3) table/spindle, with one rotation each on the table and the spindle. The generalized kinematics of 5-axis machine tools is constructed by modeling the machine as a set of links connected by the joints in an open chain. Assume that the i i i iO X Y Z and the 1 1 1 1i i i iO X Y Z+ + + + are coordinate systems attached to the links i and i+1, respectively and both links are part of the open chain. The relative position and orientation of the i+1th coordinate system (frame) with respect to the ith frame is expressed by a 4 4× homogeneous transformation matrix as: 1 1Trans( ) Rot( , ) i i i i i iN θ+ += ⋅H L (3.5) where 1, 3 3 1, 1 1, 1 3 4 4 Trans( ) 1 i ii x i ii yi ii i ii z L L L L + × + + + × × ⎡ ⎤⎢ ⎥⎢ ⎥= ⎢ ⎥⎢ ⎥⎢ ⎥⎣ ⎦ I 0 (3.6) is a homogeneous translation matrix, and 3 33 3 × × ∈ℜΙ is an identity matrix. 1 1, 1, 1, i i i i ii i i x ii y ii zL L L+ + + += + +L i j k is the relative distance from the origin iO to the origin 1iO + measured in the i th frame, as denoted by the superscript i. In order to simplify the representation, herein the superscript is omitted, and the relative displacement is always defined in the frame denoted by the first element of the subscript, such as: 1 1 i ii ii+ +=L L . (3.7) Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 24 In Equation (3.5), N represents the axis of rotation, iX , iY or iZ , where iθ is the relative rotation; hence, the pure rotation matrix representing the current frame rotation of i+1 with respect the frame i around iX yields: 1 0 0 0 0 cos( ) sin( ) 0 Rot( , ) 0 sin( ) cos( ) 0 0 0 0 1 i i i i i i X θ θθ θ θ ⎡ ⎤⎢ ⎥−⎢ ⎥= ⎢ ⎥⎢ ⎥⎣ ⎦ . (3.8) In order to describe the relative position and orientation of the cutting tool with re- spect to the workpiece, appropriate coordinate systems should be established. Figure 3.2 illustrates the placement of the coordinate systems for the table tilting 5-axis machine tool in Figure 3.1b. The coordinate systems for the workpiece and cutting tool are w w w wO X Y Z and t t t tO X Y Z , respectively. A base coordinate system, b b b bO X Y Z is fixed on the machine bed. The coordinate frames a a a aO X Y Z and c c c cO X Y Z are attached to the primary and secondary rotary joints, namely the A and C drives. The A drive is the primary axis, and C drive is the secondary rotary axis, placed closer to the workpiece. Similarly, the coordinate systems Tx Tx Tx TxO X Y Z , Ty Ty Ty TyO X Y Z and Tz Tz Tz TzO X Y Z are attached to the cartesian, X, Y, Z axes of the machine tool. Additional coordinate sys- tems, ,fw fw fw fw ft ft ft ftO X Y Z O X Y Z are assigned to model tilts or offsets between coordinate systems, such as fixtures placed between the workpiece and the table, any inclination between the rotary drives, or fixtures placed between the tool-holder and the spindle. Hence, by attaching coordinate systems both in the table and the spindle side, a genera- lized model can be developed to model kinematics of different machine tools that have rotary drives on the table and/or the spindle side. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 25 3.2.1 Forward Kinematics Modeling The forward kinematics describes the motion of the tool in the WCS as a function of the joint positions [ ], , , , Tx y z a c=Q . The generic kinematic model for 5-axis machines relies on expressing the relative position and orientation of the tool coordinate sys- tem ( )t t t tO X Y Z with respect to the workpiece coordinate system ( )w w w wO X Y Z . Both, tool and the workpiece frames exhibit spatial motion with respect to the same base coordinate system, b b b bO X Y Z . The representation of the workpiece coordinate system in base frame Zb Xb Yb ZTx XTx YTx ZTz XTy YTy Za Zc Zfw Zw θa θc θft θfw Xa Xc Xfw Xw LTy,a La,c La,fw Lfw,w LTz,ft Lft ,Tz ttp tto Ya Yc Yf1 Yw YTz XTz ZTz Yft Xft Zft Yt Xt Zt x y z workpiece tool Figure 3.2: Coordinate systems for modeling 5-axis machine tool kinematics. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 26 is obtained from Figure 3.2 by series multiplication of the homogeneous transformation matrices as: , , , , , , , , , 1 , , , Trans( ,0,0) Trans( ,0,0) Trans( , , ) Rot( , ) Trans( , , ) Rot( , ) Trans( , , ) Rot( , ) Trans( , , ) b w Ty a x Ty a y Ty a z a a a c x a c y a c z c c c fw x c fw y c fw z fw f fw w x fw w y fw w z H x y L L L X L L L Z L L L Y L L L θ θ θ = ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ (3.9) where a aθ = , c cθ = are rotations of the A and C drives, fwθ is the offset rotation of the workpiece fixture, and x, and y denote the relative translation of the X-Y table. In order to construct the generalized kinematic model, the tool frame must be represented with respect to the base frame similarly as: , , , , , , Trans( ,0,0) Trans( , , ) Rot( , ) Trans( , , ) b t Tz ft x Tz ft y Tz ft z ft ft ft t x ft t y ft t z H z L L L Z L L Lθ = ⋅ ⋅ ⋅ (3.10) where z is the translation of the Z axis. The tool tip position ( )tpt and the orientation ( )tot are defined in the tool frame: , , , , , , t t t t p p x p y p z t t t t o o x o y o z t t t t t t ⎫= + + ⎪⎬= + + ⎪⎭ t i j k t i j k (3.11) and transformed into the workpiece frame by using Equations (3.9) and (3.10) as; ( ) 1 1 0 1 0 o o w w t t p pb b w t −⎡ ⎤ ⎡ ⎤=⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦ t t t t H H , (3.12) which can be expanded to: Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 27 , , , , , 1 , , , , , , , , , Trans( , , ) Rot( , ) 1 0 Trans( , , ) Rot( , ) Trans( , , ) Rot( , ) Trans( , , ) Trans( ,0,0) Trans( ,0,0) Trans o w w p w fw x w fw y w fw z fw fw fwc x f wc y fwc z c c c a x c a y c a z a a aTy x aTy y a Ty z L L L Y L L L Z L L L X L L L y x θ θ θ ⎡ ⎤ = ⋅ −⎢ ⎥⎢ ⎥⎣ ⎦ ⋅ ⋅ − ⋅ ⋅ − ⋅ ⋅ − ⋅ − ⋅ t t , , , , , , , , , , ( ,0,0) Trans( , , ) Rot( , ) Trans( , , ) 1 0 o Tz ft x Tz ft y Tz ft z ft ft t t p ft t x ft t y ft t z z L L L Z L L L θ⋅ ⋅ ⎡ ⎤⋅ ⎢ ⎥⎢ ⎥⎣ ⎦ t t . (3.13) Equation (3.13) is called the form-shaping function [66] of the 5-axis machine. Depend- ing on the configuration, the cutting tool is generally aligned with one of the principle axis of the tool frame, namely , or t t tX Y Z . For the vertical 5-axis machining center illustrated in Figure 3.1b, the cutting tool is aligned with tZ , and its orientation is defined by the directional cosines in the augmented vector: [ ]0 0 1 0 Tto =t . (3.14) Utilizing the rotary drives of the 5-axis machine, the cutter orientation in the tool frame is aligned with the tool axis orientation in the workpiece coordinates by Equation (3.13): , , , , , , , , , , , , 0 0 Trans( , , ) Rot( , ) Trans( , , ) Rot( , ) Trans( , , ) Rot( , ) Trans( , , ) Trans( ,0,0) Trans( o o w t w b b t w fw x w fw y w fw z fw fw fw c x fw c y fw c z c c c a x c a y c a z a a a Ty x a Ty y a Ty z L L L Y L L L Z L L L X L L L y θ θ θ ⎡ ⎤ ⎡ ⎤=⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦ = ⋅ − ⋅ ⋅ − ⋅ ⋅ − ⋅ ⋅ − ⋅ − t t H H , , , , , , ,0,0) Trans( ,0,0) Trans( , , ) Rot( , ) Trans( , , ) 0 o Tz ft x Tz ft y Tz ft z ft ft t ft t x ft t y ft t z x z L L L Z L L L θ ⋅ ⋅ ⋅ ⎡ ⎤⋅ ⎢ ⎥⎢ ⎥⎣ ⎦ t (3.15) Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 28 and excluding the translation components results in =Rot( , ) Rot( , ) Rot( , ) Rot( , ) 0 0 o o w t fw fw c c a a ft ftY Z X Zθ θ θ θ ⎡ ⎤ ⎡ ⎤− ⋅ − ⋅ − ⋅⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦ t t . (3.16) Assuming that fw fw fw fwO X Y Z and ft ft ft ftO X Y Z are aligned with the previous frames, 0fw ftθ θ= = , Equation (3.16) is simplified to , , , cos( ) sin( ) 0 0 1 0 0 0 0 sin( ) cos( ) 0 0 0 cos( ) sin( ) 0 0 0 0 1 0 0 sin( ) cos( ) 0 1 0 0 0 1 0 0 0 1 00 o o o w i c c w j c c a a w a ak θ θ θ θ θ θ θ θ ⎡ ⎤ ⎡ ⎤ ⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥−⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥=⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥−⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎣ ⎦ ⎣ ⎦ ⎣ ⎦⎣ ⎦ t t t , (3.17) which yields , , , sin( )sin( ) cos( )sin( ) cos( ) 0 00 o o o w ii c a w jj c a w k ak O O O θ θ θ θ θ ⎡ ⎤⎡ ⎤ ⎡ ⎤⎢ ⎥⎢ ⎥ ⎢ ⎥−⎢ ⎥⎢ ⎥ ⎢ ⎥≡ =⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎣ ⎦ ⎣ ⎦⎣ ⎦ t t t . (3.18) Equation (3.18) proves that the tool orientation is controlled purely by the rotary drives of the 5-axis machine tool. Similarly, the tool tip position ( )tpt in the tool frame can be represented by the augmented vector as: , , , 1 Tt t t t p p x p y p z⎡ ⎤= ⎣ ⎦t t t t , (3.19) and transformed into the workpiece frame using the form-shaping function (Equation (3.13)) as: Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 29 , , , , , , , , , , , , , , , Trans( , , ) Rot( ,0) 1 1 Trans( , , ) Rot( , ) Trans( , , ) Rot( , ) Trans( , , ) Trans( w x p x w y p y w fw x w fw y w fw z fww z p z fw c x fw c y fw c z c c c a x c a y c a z a a aTy x aTy y aTy z P t P t L L L Y P t L L L Z L L L X L L L θ θ ⎡ ⎤⎡ ⎤ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ≡ = ⋅⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎣ ⎦ ⎣ ⎦ ⋅ ⋅ − ⋅ ⋅ − ⋅ ⋅ − , , , , , , ,0,0) Trans( ,0,0) Trans( ,0,0) Trans( , , ) Rot( ,0) Trans( , , ) 1 Tz ft x Tz ft y Tz ft z ft t p ft t x ft t y ft t z y x z L L L Z L L L ⋅ − ⋅ ⋅ ⋅ ⎡ ⎤⋅ ⎢ ⎥⎣ ⎦ t . (3.20) Equation (3.20) reveals that the cutter location in WCS is actually controlled by all the five physical drives of a 5-axis machine tool. 3.2.2 Inverse Kinematics Solution For a desired tool axis orientation ( ), ,i j kO O O , the rotary drive positions, aθ and cθ are solved from Equation (3.18), arccos( ) arctan2( , ) a k c i j a O c O O θ θ = = ⎫⎪⎬= = − ⎪⎭ (3.21) where arctan2( , )i jO O returns angles in the range of cπ θ π− ≤ ≤ by examining the sign of both iO and jO . Equation (3.20) is modified to derive the three Cartesian axis positions ( , , )x y z , Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 30 , , , 1 , 1 , 1 , , , , , , , , Trans( , , ) Rot( ,0) 1 Trans( , , ) Rot( , ) Trans( , , ) Rot( , ) Trans( , , ) ( 1, 1,1,1) Trans( , , ) Trans( , w fw x w fw y w fw z fw f c x f c y f c z c c c a x c a y c a z a a aTx x aTx y a Tx z Tz ft x Tz f L L L Y L L L Z L L L X L L L diag x y z L L θ θ ⎡ ⎤ = ⋅⎢ ⎥⎣ ⎦ ⋅ ⋅ − ⋅ ⋅ − ⋅ ⋅ − − ⋅ ⋅ P , , , , , , ) Rot( ,0) Trans( , , ) 1 t y Tz ft z t p ft ft t x ft t y ft t z L Z L L L ⎡ ⎤⋅ ⋅ ⎢ ⎥⎣ ⎦ t (3.22) where the translations are stacked together, and ( 1, 1,1,1)diag − − is a 4 4× diagonal matrix. Defining [ ], , Tx y z=K and assuming that 0fw ftθ θ= = for this kinematic confi- guration (see Figure 3.1b), Equation (3.22) is rewritten as: , , , , , , , , , , , , , , , , , , Trans( , , ) Trans( , , ) 1 Rot( , ) Trans( , , ) Rot( , ) Trans( , , ) ( 1, 1,1,1) Trans( , , ) 1 w fw x w fw y w fw z fwc x fwc y fwc z c c ca x ca y ca z a a a Ty x a Ty y a Ty z Tz ft x Tz ft y Tz ft z L L L L L L Z L L L X L L L diag L L L θ θ ⎡ ⎤ = ⋅⎢ ⎥⎣ ⎦ ⋅ − ⋅ ⋅ − ⋅ ⋅ − − ⎡ ⎤⋅ +⎢ ⎥⎣ ⎦ P K , , ,Trans( , , ) 1 t p ft t x ft t y ft t zL L L ⎧ ⎫⎡ ⎤⎪ ⎪⋅⎨ ⎬⎢ ⎥⎪ ⎪⎣ ⎦⎩ ⎭ t , (3.23) and the Cartesian axis position parameters ( )K are obtained as: , , , , , , , , , 1 , , , , , , {Trans( , , ) 1 Trans( , , ) Rot( , ) Trans( , , ) Rot( , ) Trans( , , ) ( 1, 1,1,1)} 1 Trans( , , ) w fw x w fw y w fw z fwc x fwc y fwc z c c ca x ca y ca z a a aTy x aTy y aTy z Tz ft x Tz ft y Tz ft z L L L L L L Z L L L X L L L diag L L L θ θ − ⎡ ⎤ =⎢ ⎥⎣ ⎦ ⋅ ⋅ − ⋅ ⋅ − ⎡ ⎤⋅ ⋅ − − ⎢ ⎥⎣ ⎦ − K P , , ,Trans( , , ) 1 t p ft t x ft t y ft t zL L L ⎡ ⎤⋅ ⎢ ⎥⎣ ⎦ t . (3.24) Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 31 Hence, for a given tool motion data T x y z x y zP P P O O O⎡ ⎤= ⎣ ⎦R , Equation (3.24) combined with the rotary drive angle solutions from Equation (3.21) represent the analytical inverse kinematic solution of the table-tilting 5-axis machine tool. 3.3 Minimum Time Feed optimization Problem In order to continuously interpolate the intermediate cutter motion, the sequence of pN discrete tool commands ( )R are fitted to a cubic, quintic, NURB or B-Spline path [67] with a total length of ∑− = Σ = 1 1 pN k kSS as shown in Figure 3.1a. The exact path length, kS , between the adjacent knots is approximated by the arc length of the parameterized splined tool-path, [68]. Tool position is expressed as a function of path arc length, s , as: Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 32 ( ) [ ( ), ( ), ( ), ( ), ( ), ( )]Tx y z i j ks P s P s P s O s O s O s=R , (3.25) and transformed into axis drive positions by the inverse kinematics model presented in Section 3.2.2 (Equations (3.21) and (3.24)) as, ( ) [ ( ), ( ), ( ), ( ), ( )]Ts x s y s z s a s c s=Q . (3.26) The velocity 5 1( )×Q , acceleration 5 1( )×Q and jerk 5 1( )×Q profiles of each drive are derived as, Figure 3.3: Tangential feed, acceleration and jerk profiles. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 33 2 3 ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) s s s s s s s s s s s s s s s s s s s s s s s ⎫= ⎪= + ⎬⎪= + + ⎭ s ss s sss ss s Q Q Q Q Q Q Q 3Q Q           , (3.27) where, 2 3 2 3 d ( ) d ( ) d ( )( ) , ( ) , ( ) d d d s s ss s s s s s ⎫= = = ⎬⎭s ss sss Q Q QQ Q Q (3.28) are the partial derivatives of axis positions with respect to the arc displacement ( )s computed analytically or numerically [69]. The derivatives, d / ds s t= , 2 2d / ds s t= and 3 3d / ds s t= are the tangential feed, acceleration and jerk profiles, 22 2 2 ( ) d ( ) d d ( )( ) ( ) d dt d d ( ) d d ( ) d ( )( ) ( ) ( ) d dt d d s s s s s s s ss s s s s s s s s s s s ss s s s s s s s s = ⎫⎪⎪= = ⎪⎬⎪⎛ ⎞ ⎪= = + ⎜ ⎟ ⎪⎝ ⎠ ⎭           , (3.29) defined along the tool-path as shown in Figure 3.3. Equation (3.27) implies that for a given tool-path, the axis kinematic profiles, velocity ( )Q , acceleration ( )Q and jerk ( )Q are controlled by the tangential kinematics of the trajectory ( ), ,s s s   . The objective of feed optimization is to minimize 5-axis machining cycle time by scheduling the feed, ( )s s s=  , along the path displacement in such a way that the con- straints are not violated. The instantaneous feed along the path is /s ds dt= ; equivalently, the differential time ( )dt to travel along infinitesimal path segment ( )ds is d d /t s s=  . The feed optimization problem is defined as the minimization of total travel Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 34 time, TΣ , along the entire path, SΣ , while respecting a set of physical limits, C , of the machine and the process: { } { } 0 dmin( ) min subject to , , Axis Limits for 0 s S ( ) S s sT s s s C s s Σ Σ Σ= ≤ ≤ ≤∫    . (3.30) An arc length parameterized cubic B-spline is proposed in this work to represent a jerk continuous )( 3C feed profile as a function of path position as: ( ) ( ) ( ) ( )10,3 0 1,3 1 1,3 1 ,3 0 ( ) n n n i i i s s N s f N s f N s f N s f − − − = = ⋅ + ⋅ + + ⋅ = ∑ … (3.31) where , ( 0,1... 1)if i n= − are the n number of feed control points fixed along the path length 0 1 1[ , ,..., ] p p p p ns s s −=s , as shown in Figure 3.3. The nominal feed speeds are assigned to the control points are interpolated by a set of cubic blending functions ( )sNi 3, to construct the B-spline feed profile along the tool-path. By increasing the nominal feedrate assigned to a control point, if , the path velocity in the vicinity of p is is increased, and the entire feed profile can be optimized by modulating all the control points. Blending functions are computed discretely through evaluation of the known recurrence formula [70], ( ) ( ) ( ) ( ) 1 ,0 1 , , 1 1, 1 1 2 1 1 0 . , 1,2,3. i i i i i i k i k i k i i i i if u u u N u otherwise k u u u uN u N u N u u u u u + + − + − + + + ⎫≤ ≤⎧= ⎨ ⎪⎪⎩ =⎬− − ⎪= + ⎪− − ⎭ , (3.32) where [ ]0,1u∈ is the normalized path parameter computed as, Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 35 1u s SΣ ≅ ⋅ . (3.33) The vector u parameterizes the feed profile, computed through uniformly distributing n- 2 parameters along the path as: 1 2 n-2 points 1[0 ... 1] [0 2 ... ] , where 3 Su u s s S s S n Σ Σ Σ = = ∆ ∆ ⋅ ∆ = −u  . (3.34) In order to force the feed profile to pass through the initial and the final feed constraints 00, inits f f= = and 1, n finals S f fΣ −= = , the first and last entries are repeated yielding the augmented parameter vector: N N1 4 4 [0...0 ... 1...1]i iu u +=u . (3.35) The 1st and 2nd derivatives of the feed profile are obtained from Equation (3.31) and (3.33) as, ( ) ( ) ( ) ( ) ( ) ( ) ( ) 1 0,3 1,3 ,3 0 1 0 2 2 22 1 0,3 1,3 ,3 0 12 2 2 2 0 d d dd ( ) d d d d d dd ( )( ) d d d d n n i n i i n n i n i i N u N u N us u s u f f f u u u du N u N u N us us u f f f u u u u −− − = −− − = ⎫⎛ ⎞= = ⋅ + + ⋅ = ⎪⎜ ⎟⎝ ⎠ ⎪⎬⎛ ⎞⎪= = ⋅ + + ⋅ = ⎜ ⎟⎪⎝ ⎠⎭ ∑ ∑  …  … (3.36) and the derivatives of the blending functions are evaluated through repeated differentia- tion, ( ) ( ) ( )(k) (k-1) (k-1) (k-1) (k-1)3 2 1 2(k) 3 4 1 d d d d d 3 , 1,2 d i, i, i , i i i i N u N u u N u u k u u u u u + + + + ⎛ ⎞= − =⎜ ⎟− −⎝ ⎠ . (3.37) Hence, the position dependent tangential, acceleration and jerk profiles are obtained by de-normalizing Equation (3.36) using Equation (3.33), Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 36 22 2 2 2 1 d ( )( ) ( ) 1 d ( ) 1 d ( )( ) ( ) ( ) d d u s S u s Su s S s us s s u S du s u s us s s u s u S u S u Σ ΣΣ Σ = Σ Σ == ⎫⎛ ⎞⎛ ⎞= ⎪⎜ ⎟⎜ ⎟⎝ ⎠⎝ ⎠ ⎪⎪⎬⎛ ⎞⎛ ⎞ ⎛ ⎞⎛ ⎞ ⎪⎛ ⎞= ⋅ + ⋅⎜ ⎟⎜ ⎟ ⎜ ⎟⎜ ⎟ ⎜ ⎟ ⎪⎜ ⎟⎝ ⎠⎝ ⎠⎝ ⎠ ⎝ ⎠ ⎪⎝ ⎠ ⎭       . (3.38) 3.4 Kinematic Constraints The user commands the tool to maintain the feed along the tool-path in order to pre- serve the allowed chip load on the tool, i.e., feeds ≤ . However, as the curvature of the path changes, the commanded path velocity, the feed, may violate the speed, acceleration and jerk limits of active drives on the 5-axis machine tool. The feed must be constrained as follows to avoid the violation of axis limits. The following subsection addresses those constraints. 3.4.1 Velocity Constraints The velocities of the five drives as given in Equation (3.27), must not exceed their saturation limits, ,max ,max ,max ,max ,max T x y z a cv v v v v⎡ ⎤⎣ ⎦ , for a 5-axis machine tool traveling along a curved tool-path ( Σ≤≤ Ss0 ), , ax ,max ,max ,max ,max ( )1 1 ( )1 1 ( ) ( )1 1 ( )1 1 ( )1 1 s x m s y s z s a s c s x s v y s v z s v s s a s v c s v ⎡ ⎤⎡ ⎤ ⎡ ⎤⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥− ≤ ≤⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦⎣ ⎦ Q   . (3.39) The path velocity is positive 0)( ≥ss , and using the normalized partial derivatives, Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 37 5 1 ,max ,max ,max ,max ,max ( ) ( ) ( ) ( ) ( )( ) , , , , T v s s s s s s x y z a c x s y s z s a s c ss v v v v v× ⎡ ⎤= ⎢ ⎥⎢ ⎥⎣ ⎦ Q (3.40) the axis velocities are constrained, to satisfy ( ) ( ) 1 0, 1 2 5vs,iQ s s s i , ...⋅ − ≤ = . (3.41) The absolute values of the partial derivatives of the axis positions ( ), ( )vs iQ s are ex- pressed as scalar functions of the path length ( )s ; thus they are evaluated from the specified path geometry as shown in Figure 3.1a. The minimum of allowable feeds ( )Vs , computed from Equation (3.41) satisfies the velocity ( ),max ,max ,max ,max ,max Tx y z a cv v v v v⎡ ⎤⎣ ⎦ constraints of all five drives: { },max 1,2...5 ,max 1,2...5 1 1 1( ) min ( ) max ( )V vv vi ss,i s,ii s s QQ s Q s= = ⎧ ⎫⎪ ⎪= = =⎨ ⎬⎪ ⎪⎩ ⎭  . (3.42) It is evident from Equations (3.39) and (3.42) that the drive for which has the maximum first derivative of its displacement is normalized with its velocity limit, determines the allowable feed along the path ( )s . Equation (3.42) may also be expressed in a single velocity inequality as: ,max( ) ( ) 1 0 v V sC s Q s s= − ≤ . (3.43) 3.4.2 Acceleration Constraints Either the torque or the equivalent acceleration at each drive must be constrained to avoid actuator saturation limits. If the limits of five drive accelerations are specified as Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 38 ,max ,max ,max ,max ,max T x y z a ca a a a a⎡ ⎤⎣ ⎦ , the feed ( )s must be constrained according to Equation (3.27) as, ,max ,max ,max ,max 2 ,max ,max ,max ,max ,max ,max ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) x xss s y yss s z zss s a ass s c css s ss s a ax s x s a ay s y s a az s s s z s s s a aa s a s a ac s c s ⎡ ⎤ ⎡⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢⎢ ⎥ ⎢ ⎥− ≤ + ≤⎢ ⎥ ⎢⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦⎣ ⎦ ⎣ Q Q     ⎤⎥⎥⎥⎥⎢ ⎥⎢ ⎥⎦ (3.44) for Σ≤≤ Ss0 . By normalizing the partial derivatives of the axis positions with the individual axis acceleration limits: ,max ,max ,max ,max ,max ,max ,max ,max ,max ,max ( ) ( ) ( ) ( ) ( )( ) , , , , ( ) ( ) ( ) ( ) ( )( ) , , , , T a s s s s s s x y z a c T a ss ss ss ss ss ss x y z a c x s y s z s a s c ss a a a a a x s y s z s a s c ss a a a a a ⎫⎡ ⎤ ⎪= ⎢ ⎥ ⎪⎢ ⎥⎣ ⎦ ⎪⎬⎡ ⎤ ⎪= ⎢ ⎥ ⎪⎢ ⎥ ⎪⎣ ⎦ ⎭ Q Q , (3.45) Then, equation (3.44) is reduced to the following form: 2( ) ( ) ( ) ( ) 1, 1 2 5a ass,i s,iQ s s s Q s s s i , ...+ ≤ =  . (3.46) Using the triangular ( )BABA +≤+ and multiplicative ( )BABA ⋅≤⋅ properties of the absolute value, the acceleration inequalities may be expanded as: 2, ,( ) ( ) ( ) ( ) 1, 1,2...5 a a ss i s iQ s s s Q s s s i+ ≤ =  . (3.47) Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 39 Note that Equation (3.47) represents a linear relationship between square of the feed )( 2s and the tangential acceleration )(s for each axis [ ]( )( ) ( ), ( ), ( ), ( ), ( ) Ts x s y s z s a s c s=Q , as shown in Figure 3.4. When the machine starts to move from zero speed, it can use the maximum allowable tangential accelera- tion to reach the desired feed, which corresponds to ( )s axis. When the machine travels at the constant feed ( )s , the required tangential acceleration is zero. The acceleration limits of the drives determine the allowable feed and tangential acceleration. While one axis may limit the feed, another may constrain the acceleration depending on their partial derivatives of the path. A pair of feed ( )s and tangential acceleration ( )s must be identi- fied iteratively during feed optimization from the five sets of inequalities to avoid Figure 3.4: Acceleration feasible operating region. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 40 violating the acceleration limits of the drives which leads to the acceleration based feasible operating region of the machine at each path position as illustrated in Figure 3.4. However, during an iterative optimization scheme the solution of feasible feed and acceleration values may require unaffordable computational costs for the real time trajectory generation algorithm. The feasible solution region for the feed and tangential acceleration may be bounded by the acceleration limits of only a few drives as illustrated in Figure 3.4, and evaluated from: ,max 1,2...5 ,max, ,1,2...5 ,max 1,2...5 ,max, ,1,2...5 1 1 1( ) 0 ( ) min ( )( ) max ( ) 1 1 1( ) 0 ( ) min ( )( ) max ( ) A aa ai ss i s ii A aa ai ssss i ss ii s s s s Q sQ s Q s s s s s Q sQ s Q s = = = = ⎫⎧ ⎫⎧ ⎫⎪ ⎪ ⎪ ⎪ ⎪= ⇒ = = =⎨ ⎬ ⎨ ⎬ ⎪⎪ ⎪ ⎪ ⎪⎩ ⎭ ⎩ ⎭ ⎪⎬⎧ ⎫⎧ ⎫ ⎪⎪ ⎪ ⎪ ⎪= ⇒ = = =⎨ ⎬ ⎨ ⎬⎪ ⎪ ⎪ ⎪⎩ ⎭ ⎩ ⎭ ⎭     ⎪⎪ . (3.48) By substituting ,max ,max, a a s ssQ Q , Equation (3.47) can be reduced to a simple linear relation- ship defining a triangular, feasible feed and tangential acceleration region bound by ,max ,max,A As s  , 2,max ,max( ) ( ) ( ) ( ) ( ) 1 0 a a A s ssC s Q s s s Q s s s= + − ≤  (3.49) where the acceleration limited maximum reachable feedrate at a point of the path *s s= is found by setting the acceleration *( ) 0s s = in Equation (3.49), *,max * ,max 1( ) ( )A ass s s Q s = . (3.50) Since *,max ( )As s is calculated for zero path acceleration ( )*( ) 0s s = at that particular point, the tangential feed cannot be altered as the machine does not have available Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 41 acceleration due to the saturation limits for braking to a lower velocity or increasing to a higher speed as the tool continues to travel. Therefore, set of feasible feed and accelera- tion along the entire length of the tool-path must be identified in an iterative fashion, which is solved utilizing the optimization algorithm presented in Section 3.5. Equation (3.49), on the other hand, reduces the computational cost of the iterative feed optimiza- tion procedure by cutting down the number of constraint evaluations from five (Equation (3.44)) to one, at the expense of losing some part of the feasible feed-acceleration region because of the triangular inequality (Equation (3.47)) and the max operations (Equation (3.48)). 3.4.3 Jerk Constraints Limiting the axis jerk values smooths the position commands, hence the high fre- quency content, which leads to a better tracking performance along the curved path. If the imposed axis jerk limits are ,max ,max ,max ,max ,max T x y z a cJ J J J J⎡ ⎤⎣ ⎦ , the following norma- lized jerk constraints are defined from Equation (3.27): ,maxmax ,maxmax 3 ,maxmax ,maxmax ,maxmax ( )( )1 ( )( )1 ( )( ) ( ) 31 ( )( )1 ( )( )1 ss xsss x ss ysss y ss zsss z ss asss a ss csss c j j sss ss x s Jx s J y s Jy s J z s Jz s J s s a s Ja s J c s Jc s J ⎡ ⎤⎡ ⎤⎡ ⎤ ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥− ≤ + ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥⎣ ⎦ ⎣ ⎦ ⎣ ⎦ Q Q    ,max ,max ,max ,max ,max ( ) 1 ( ) 1 ( )( ) ( ) ( ) 1 ( ) 1 ( ) 1 s x s y s z s a s c j s x s J y s J x s Js s s s s s x s J x s J ⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥+ ≤⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦⎣ ⎦ Q     , (3.51) for Σ≤≤ Ss0 . Alternatively, Equation (3.51) can be expressed by, 3, , ,( ) ( ) 3 ( ) ( ) ( ) ( ) ( ) 1, 1,2...5 j j j sss i ss i s iQ s s s Q s s s s s Q s s s i+ + ≤ =    , (3.52) Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 42 and re-written as: 3, , ,( ) ( ) 3 ( ) ( ) ( ) ( ) ( ) 1, 1,2...5 j j j sss i ss i s iQ s s s Q s s s s s Q s s s i+ + ≤ =    , (3.53) where, ( )jsss,iQ s , ( ) j ss,iQ s and , ( ) j s iQ s are tool-path geometry dependent terms. Equation (3.53) generates a three dimensional operating surface with respect to the jerk limits of the drives as shown in Figure 3.5. Similar to the acceleration constraints presented in Section 3.4.2, the axis jerk limits transform into upper bounds on the tangential feed, s . By setting both the tangential acceleration ( )s and the jerk ( )s to zero at a path position, *s s= , the jerk limited maximum path velocity, ,max ( )Js s , is evaluated from Equation (3.53) as, * * * 3,max *1,2...5 , ( ) ( ) 0 ( ) min ( )J ji sss i s s s s s s Q s= ⎧ ⎫⎪ ⎪= = ⇒ = ⎨ ⎬⎪ ⎪⎩ ⎭ 1   . (3.54) * ,max ( )Js s represents an upper limit on the path velocity with respect to the jerk limits of the drives and cannot be altered, since it utilizes the entire acceleration and jerk capacity of the drives. Therefore, the feasible set of feed, tangential acceleration and jerk identi- fied in the operating region is then iteratively found during feed optimization by evaluating Equation (3.53), which is computationally complex and too costly to be implemented in a CNC system. Similar to the methodology presented for the accelera- tion limits, the maximum values of the jerk normalized partial derivatives are selected by considering all five drives: Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 43 { } { } { } ,max ,1,2...5 ,max ,1,2...5 ,max ,1,2...5 ( ) max ( ) ( ) max ( ) ( ) max ( ) j j sss sss ii j j ss ss ii j j s s ii Q s Q s Q s Q s Q s Q s = = = ⎫= ⎪⎪= ⎬⎪⎪= ⎭ (3.55) and substituted in Equation (3.53) to find the constraint equation: 3.max ,max ,max( ) ( ) ( ) 3 ( ) ( ) ( ) ( ) ( ) 1 0 j j j J sss ss sC s Q s s s Q s s s s s Q s s s= + ⋅ + − ≤   . (3.56) While Equation (3.53) has to be evaluated for all five axes, Equation (3.56) is a single set, which considers a combination of axes whose jerk limits affect the feasible operating region of the machine tool. The tangential jerk, acceleration and speed are still coupled; hence, they need to be identified iteratively during the optimization process from Equa- s s s X Axis Jerk Limits Y Axis Jerk Limits Modified Feasible Operating Region Global Jerk Constraint Figure 3.5: Jerk feasible operating region. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 44 tion (3.56). The identified feed )(ssJ , tangential acceleration )(ssJ and tangential jerk )(ssJ satisfy the jerk limits on all five drives. 3.5 Solution to Minimum Time Feed Optimization Problem It is proposed that the feed is scheduled along the tool-path without violating the ve- locity, acceleration and jerk limits of all five drives. The axis constraints given in Equations (3.43), (3.49), (3.56) are stacked into: 5 1 ( ) ( ) ( ) ( ) ( ) ( ) T F F V A Js C s C s C s C s C s − + ×⎡ ⎤= ⎣ ⎦C , (3.57) where ( ) ( ) 0FC s s s − = − ≤ is the lower limit of the velocity to ensure forward motion (i.e., positive feed direction) , and ( ) ( ) 0FC s s s feed + = − ≤ is the user specified upper limit of the feed along the path. The feed profile is expressed in B-spline form and optimized to minimize the machining time by maximizing the nominal feeds at control points: ( ) ( ) 1 min max subject to 0 for 0 S i n is f i T f C s sΣ Σ = = ≤ ≤ ≤∑ (3.58) The feed control points if are the optimization variables, and the constraints are eva- luated from Equation (3.58). Since the second derivative of TΣ with respect to if is zero, solving Equation (3.58) becomes a Linear Programming problem with nonlinear con- straints yielding a Sequential Quadratic Problem (SQP) and solved using MATLAB’s Optimization Toolbox [71]. Five to seven feed control points are distributed along each segment length kS , and the 5 ~ 7 pn N= nominal feeds are modulated while the con- straint violations are checked at denser increments [ ]( )0.1...1s mm∆ ≈ along the path Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 45 length as shown in Figure 3.6. As the feed profile is altered iteratively, a feasible set of control points are obtained that does not violate the constraints and that minimize the overall cycle time. The optimization process is computationally inefficient if all feed control points are considered at once along the entire tool-path, due to the large number of optimization variables and constraints that need to be considered. Instead, the optimization algorithm is applied in a moving window approach that considers a smaller set of control points over each window as shown in Figure 3.7. The solution window is started in the forward Figure 3.6: Feed optimization along the tool-path. Displacement (s) Displacement (s) Displacement (s)Displacement (s) a) Control Point fi Displacement interval Optimized Feed Profile Feed Control Points Extended Displacement Interval b) c)d) Displacement Interval Overlap Optimized Feed Profile ∆s Constraint check points Fe ed P ro fil e (s ) Fe ed P ro fil e (s ) Fe ed P ro fil e (s ) Fe ed P ro fil e (s ) Figure 3.7: Windowing technique: a) forward windowing (start), b) forward windowing (shifted), c)forward windowing (finalized), c) forward winding (extended). Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 46 direction with an initial number of 20n ≈ control points tuned for quick convergence. Upon convergence the window is moved along the tool-path while allowing overlap of m=5 control points from the previous solution, and the next solution is sought. If the trajectory optimization in a particular window does not converge to a feasible profile, the solution is backtracked by expanding the window size. For instance, the final feedrate at the end of the tool-path is set to be zero to demand a full stop from the machine tool. In this case, the last window is expanded by borrowing control points from previous windows until a feasible braking distance is obtained. Figure 3.8: Five-axis tool-path. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 47 3.6 Simulations and Experimental Results The proposed trajectory optimization has been tested on a open-architecture con- trolled 5-axis machine tool with the tilting-rotary table kinematics studied in Section 3.2 (see Figure 3.1b). The test tool-path shown in Figure 3.8 is designed for 5-axis flank milling of jet impellers. The discrete tool-path containing 196 data points is fitted into a C3 quintic B-spline with a tool tip chord error of 35 [micron] and an orientation error constrained under 0.09 [degree]. The curvature along the tool-path changes continuous- ly, but the highest curvatures are located at the beginning (25~50 [mm]), in the middle (150~190 [mm]) and towards the end (250~290 [mm]) where the tool has to move around the trailing and leading edges of the blade. The corresponding axis positions along the 300 [mm] long tool-path are shown in Figure 3.9. The programmed feed for the tool-path is 150 [mm/sec] in the CL file. The experimental 5-axis machine tool setup has been introduced in Chapter 1. The axis kinematic limits are given in Table 3.1, and the drive system is controlled by robust sliding mode controllers with feed-forward friction compensation [35]. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 48 Table 3.1: Drive limits used in the feed optimization. Drives Velocity [unit/sec] Acceleration [unit/sec2] Jerk[unit/sec3] X [mm] 250 2000 3000 Y [mm] 250 2000 3000 Z [mm] 250 2000 3000 A [deg] 166 2500 250 C [deg] 290 2000 200 The path velocities constrained by the drive velocity, acceleration and jerk limits are evaluated from Equations (3.42), (3.50) and (3.54), and shown in Figure 3.10. Depend- ing on the curvature of the path, different drives limit the feed along the tool-path. For example, the X axis limits the maximum feed at 58=s [mm/sec] at path position 34=s [mm], and the Y and the C drives limit the feed at path location 267=s [mm/sec]. Similarly, the acceleration and jerk limits of the drives constrain the feed at different locations of the path. The envelope, which passes through minimum feeds 0 50 100 150 200 250 300 -100 -50 0 50 100 Displacement[mm] [m m ] 0 50 100 150 200 250 300 -6 -5 -4 -3 -2 -1 0 Displacement[mm] [r ad ] x Axis y Axis z Axis a Axis c Axis Cartesian axis displacement,{x(s), y(s), z(s)} Rotary axis displacement, {a(s), c(s)} Figure 3.9: Axis reference displacement. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 49 limited by the velocity (Figure 3.10a), acceleration (Figure 3.10b) and jerk (Figure 3.10c) is evaluated along the path. The envelope curve can be taken as a guide for selecting feasible feeds along the toolpath. However, it is not simply possible to accept the feedrate defined by the slowest envelope curve as the optimal velocity profile. Since a change in the velocity between any two points along the envelope curve may violate the acceleration or the jerk limits of the drives, the optimal feed must be evaluated by simultaneously optimizing velocity, acceleration and jerk limits of all five drives as given in Equations (3.43), (3.49) and (3.56). The identified feed is limited to 150 [mm/sec], set by the process planner, and the resulting feed schedule along the path is given in Figure 3.11. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 50 It can be seen that, if the machine travels at 150 [mm/sec], the drive limits are vi- olated at path regions (0-60 [mm], 140-165 [mm], 240-300 [mm]), where the impeller blade curvature is high. If the feed is kept constant at 15 [mm/sec], none of the drive limits are violated, but the machining time becomes [sec]18.20=ΣT . When the feed is optimally scheduled by fitting a B-spline to n=276 equally spread points along the tool- path, the machining time drops to [sec]664.4=ΣT , resulting in an 80% increase in 0 50 100 150 200 250 3000 200 400 600 800 1000 1200 Displacement[mm] [m m /s ec ] 0 50 100 150 200 250 3000 200 400 600 800 1000 1200 0 50 100 150 200 250 3000 200 400 600 800 1000 1200 Displacement [mm] a) b) c) Displacement [mm] X Axis Limit Y Axis Limit Z Axis Limit A Axis Limit C Axis Limit [m m /s ec ] [m m /s ec ] Figure 3.10: Axis limits reflected on the tangential feed: a) feedrate limited by the axis velocity constraints, b) feedrate limited by the axis acceleration constraints, c) feedrate limited by the axis jerk constraints. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 51 productivity. The algorithm was solved in the MATLAB® environment in 126 [sec] on a P4 computer with 2.3 Ghz CPU speed. As shown in Figure 3.12, the scheduled feed has been tested on the experimental 5-axis machine, and the resulting velocity, accelera- tion and jerk on each drive did not violate their limits. Moreover, depending on the complexity of 5-axis tool-paths different computational times have been observed for solving the optimization problem. However, the results did not violate the limits im- posed in the algorithm. In addition to the reduction in machining time, the proposed optimization improves the tracking accuracy of the machine when the drives are not saturated. The experimen- tally measured tracking errors are shown in Figure 3.13 for three cases. When the feed is set to the most conservative speed (15 [mm/sec]), the drive limits are never violated, and the tracking accuracy is best due to the slow speed, but the machining time be- comes [sec]18.20=ΣT . If the user imposed feed (150 [mm/sec]) had used, the drives would have become excessively saturated and the machine would become uncontrolla- ble. A relatively conservative feed of 40=f [mm/sec] was tested, leading to the violation of drive limits (Figure 3.13), and the tracking errors were amplified at the 0 50 100 200 250 3000 50 100 150 200 250 Velocity Constraints Acceleration Constraints Jerk Constraints Displacement(s) [mm] Fe ed ra te (s ) [ m m /s ec ] Constant Feed Profile Optimized Feed Profile Feed Limit Figure 3.11: Optimized and un-optimized feed profiles. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 52 saturation regions along the path. However, when the feed is optimized, the machine was able to reach up to the desired feed of 150 [mm/sec] at most locations, and the drive limits were never violated. The machine tool control was kept in a linear control region and the tracking accuracy became linearly proportional to the velocities and inversely proportional to the bandwidth of the drives. 0 4 8 16 20 -5 5 0 4 8 16 20 -50 50 -300 300 A Axis Limit C Axis Limit A Axis C Axis Rotary axis acceleration profiles Rotary axis jerk profiles A Axis Limit C Axis Limit 0 4 8 16 20 -200 0 0 200 0 4 8 88 16 20 2x103 -2x103 3x104 -3x104 0 X Axis Y Axis Cartesian axis velocity profiles (Optimized profiles (solid), un-optimized (dashed) Cartesian axis acceleration profiles Z Axis Cartesian axis jerk profiles 0 4 16 20 0 4 16 20Time [sec]Time [sec] Time [sec]Time [sec] Time [sec]Time [sec] [m m /s ec ] [m m /s ec 2 ] [m m /s ec 3 ] [d eg /s ec 3 ] [d eg /s ec 2 ] [d eg /s ec ] X Axis Y Axis Z Axis X Axis Y Axis Z Axis X Axis Y Axis Z Axis A Axis C Axis A Axis C Axis A Axis C Axis Rotary axis velocity profiles (Optimized profiles (solid), un-optimized (dashed) Limit Limit Limit Limit Figure 3.12: Optimized and un-optimized axis profiles. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 53 3.7 Conclusion The productivity of 5-axis machining of dies, molds and aerospace parts with sculp- tured surfaces is hindered by the metal cutting process, the machine tool’s structural dynamics, thermal and volumetric accuracy of the machine, and kinematics and control of the drives. The algorithm developed in this chapter has significantly reduced the 0 5 15 20 0 0 5 15 20 0 0 5 15 20 0 0 5 15 20 0 5 15 20 Tracking errors of X axis Tracking errors of Y axis Tracking errors of Z axis Tracking errors of A axis Tracking errors of C axis 0.5 -0.5 0.5 -0.5 0.5 -0.5 0.5 -0.5 0.5 -0.5 Time [sec] Time [sec] Time [sec] Time [sec] Time [sec] F=150[mm/sec](Optimized) [m m ] [m ill i- de gr ee ] [m ill i- de gr ee ] [m m ] [m m ] F=15[mm/sec](Const.) F=40[mm/sec](Const.) Actuator Saturation Actuator Saturation Actuator Saturation Figure 3.13: Tracking errors of the axes. Chapter 3: Feed Optimization for 5-Axis CNC Machine Tools with Drive Constraints 54 machining time by optimizing the feed along the tool-path while respecting velocity, acceleration and jerk limits of all five drives and reported in [15]. In addition, by avoid- ing violation of the saturation limits of the drives, the control law can operate in a linear region, improving the accuracy and performance of the CNC system. Moreover, the tangential feed along the cutter axis varies when the axial depth of cut is long during 5- axis flank milling, leading to conflicting velocity, acceleration and jerk demands at each drive. The proposed algorithm can also be repeated along the cutter axis by considering the tangential feed at that position, and a feed which satisfies the drive constraints can be selected. 55 4. Modeling of Contouring Error for 5-Axis Machine Tools 4.1 Introduction Each drive on a 5-axis machine tool is subject to unavoidable tracking errors, , , , ,x y z a ce e e e e , between the commanded axis positions, , , , ,ref ref ref ref refx y z a c , and the actual axis positions, , , , ,x y z a c , due to the servo dynamics. Axis tracking errors are kinematically transformed onto the tool-path, and the normal deviation of the actual tool trajectory from the commanded tool-path geometry is denoted by the contour errors. The contour errors control the dimensional accuracy of the finished part during 5-axis high speed machining. Due to increase in the demanded feedrate, higher axis tracking errors are observed on the drives. This in return, may induce more contour errors that are left on the part surface. In the current practice, dedicated, high-bandwidth, control laws are designed for in- dividual axis drives to accurately track axis reference commands. On the other hand, as a consequence of the non-Cartesian kinematics of the 5-axis machine tool and the complex tool-path geometry, small axis tracking errors may result in larger deviation of the tool from the reference tool-path relative to the workpiece. Correctly estimating the resultant contour errors in real time, and subsequently minimizing them with a robust controller, improves the 5-axis machining accuracy. This chapter presents a novel technique for definition and modeling of contour errors that occur during simultaneous 5-axis machin- ing. During 5-axis motion the tool translates and simultaneously changes its orientation (Figure 4.1a). The normal deviation of the tool tip position from the reference tip trajec- tory is modeled, and is denoted as the tool tip contour error. The tool tip contour error, Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 56 illustrated in Figure 4.1b, is crucial in avoiding both over-cutting and under-cutting during general point milling operations, and analytical models to accurately estimate tool tip contour errors are presented in Section 4.3. On the other hand, the normal deviation of the actual tool orientation form the desired orientation trajectory is vital when the depth of cut is high, and the cutting is mainly performed along the tool axis, as in flank milling of impeller blades. The orientation contour error for 5-axis machining is shown in Figure 4.1b, and models to accurately predict them are presented in Section 4.4. Both contour error estimation models are validated experimentally in Section 4.5, and the contributions of the chapter [26] are presented in Section 4.6. P = [Px Py Pz] T O = [Oi Oj Ok] T X Y Z Workpiece Coordinate System Cutter Contact Point Tool Tip Point a) b) ε Actual Tool Tip PointSurface φ Actual Tool Orientation Figure 4.1: Definition of contour errors in 5-axis machining: a) tool motion in work- piece coordinates (P-system), b) contour error definition (orientation and tool tip contour errors). Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 57 4.2 Tracking Error Dynamics on the Workpiece Contour errors left on the finished surface vary as a non-linear function of the track- ing errors of the tool (Figure 4.1a). Tracking errors of the tool, on the other hand, are governed by both the kinematics of the machine tool and the tracking errors of the axis servo drives ( , , , ,x y z a ce e e e e ). Therefore, axis tracking errors are transformed to the tool- path using kinematic model of the machine to prodict the actual tool trajectory. A widely used tilting-rotary-table type kinematic configuration has been presented in Chapter 3. The 5-axis machine has three translational motions, X, Y, Z, combined with rotations of the A and C drives. The transformation between the tool motion , , , , , T x y z i j kP P P O O O⎡ ⎤= ⎣ ⎦R in the WCS and the motion of the physical drives, [ ], , , , Tx y z a c=q , has also been presented in Chapter 3. The angular motions of the rotary axes were obtained as arccos( ) arctan2( , ) a k c i j a O c O O θ θ = = ⎫⎪⎬= = − ⎪⎭ , (4.1) and the forward kinematics yields the tool axis direction as: sin( )sin( ) cos( )sin( ) cos( ) i c a j c a k a O O O θ θ θ θ θ ⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥= −⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦ . (4.2) The rotary axis velocities ( ), ( )a t c t  are then related to the angular velocity vector, ( ) [ ( ), ( ), ( )]Ti j kt O t O t O tω =G    , of the tool by the orientation Jacobian 3 2×∈ℜoJ of the machine obtained from Equation (4.2): Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 58 N sin( )cos( ) cos( )sin( ) cos( )cos( ) sin( )sin( ) sin( ) 0 i c a c a a j c a c a c k a O O O ω θ θ θ θ θθ θ θ θ θθ ⎡ ⎤ ⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥= − − ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎣ ⎦⎢ ⎥ ⎢ ⎥−⎣ ⎦⎣ ⎦ Jo G     . (4.3) Equations (4.2) and (4.3) reveal that only the two rotary drives, A and C, control the tool orientation ( ), ,i j kO O O . Therefore, the tracking errors of the rotary drives, ,a ce e will affect the accuracy of the tool orientation. The forward kinematics define the tool tip trajectory, , , T x y zP P P⎡ ⎤= ⎣ ⎦P as a function of the axis motion, , , , , , , , , , , , , Trans( , , ) Rot( ,0) 1 1 Trans( , , ) Rot( , ) Trans( , , ) Rot( , ) Trans( , , ) Trans( ,0,0) Trans( ,0,0) x y w fw x w fw y w fw z fw z fw c x fw c y fw c z c c c a x c a y c a z a a aTy x aTy y aTy z P P L L L Y P L L L Z L L L X L L L y x θ θ ⎡ ⎤⎢ ⎥ ⎡ ⎤⎢ ⎥ = ⋅⎢ ⎥⎢ ⎥ ⎣ ⎦⎢ ⎥⎣ ⎦ ⋅ ⋅ − ⋅ ⋅ − ⋅ ⋅ − ⋅ − ⋅ P , , , , , , Trans( ,0,0) Trans( , , ) Rot( ,0) Trans( , , ) 1 Tz fp x Tz fp y Tz fp z ft t p fp t x fp t y fp t z z L L L Z L L L ⋅ ⋅ ⎡ ⎤⋅ ⎢ ⎥⎣ ⎦ t , (4.4) and we can simplify this expression by setting some of the x and y offsets to zero at the home position: Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 59 cos( ) sin( )cos( ) sin( )sin( ) sin( ) cos( )cos( ) cos( )sin( ) 0 sin( ) cos( ) sin( )sin( ) sin( )sin( ) 0 + cos( )sin( ) cos( )sin( ) 0 cos( ) c x c c a c a y c c a c a z a a c a c a c a c a a P x P y P z θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ −⎡ ⎤ ⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥ ⎢ ⎥= − − −⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥−⎣ ⎦ ⎣ ⎦ ⎣ ⎦ − − − , , ,os( ) 1 Tya z fpt z a ac z L L Lθ ⎡ ⎤⎡ ⎤ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥−⎣ ⎦ ⎣ ⎦ . (4.5) The inverse kinematic transformation yields the Cartesian axis positions as: , , , , cos( ) sin( ) 0 cos( )sin( ) -cos( )cos( ) sin( ) sin( )sin( ) -sin( )cos( ) cos( ) 0 sin( ) cos( ) + - c c x a c a c a y a c a c a z a ac z a ac z Tya z fpt z x P y P z P L L L L θ θ θ θ θ θ θ θ θ θ θ θ θ θ − −⎡ ⎤ ⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥ ⎢ ⎥= −⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦ ⎣ ⎦ ⎡ ⎤⎢ ⎥+ −⎢ ⎥⎢ ⎥⎣ ⎦ . (4.6) In contrast to the orientation tracking, the direct kinematic transformation, Equation(4.5) , shows that translational (X, Y, Z) as well as the motion of the rotary drives (A,C) govern the tool tip trajectory ( )[ ]Tx y zP P P P= on the workpiece surface. Thus, tracking errors of all the drives ( ), , , ,x y z a ce e e e e are reflected on the tool tip trajectory and restrain to the positioning accuracy of the tool tip. A ‘tool pose’ variable is defined as: ( ) [ ( ), ( ), ( ) ( ), ( )]Tx y zt P t P t P t a t c t= P x  , (4.7) containing the tool tip location in WCS and the two rotary drive positions. The partial derivatives of tool tip ( )x y zP P P with respect to tool pose variables ( ), , , ,x y zP P P a c are evaluated by differentiating Equation (4.5), Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 60 ( , , ), ( , , , , )iij j j PJ i x y z q x y z a c q ∂= ∈ ∈∂ . (4.8) By adding two unity elements on the diagonal for the rotational drive positions, the following square Jacobian matrix ( )5 5( )t ×∈ℜJ is formed: 11 12 13 14 15 21 22 23 24 25 31 32 33 34 35( ) 0 0 0 1 0 0 0 0 0 1 J J J J J J J J J J t J J J J J ⎡ ⎤⎢ ⎥⎢ ⎥⎢ ⎥= ⎢ ⎥⎢ ⎥⎢ ⎥⎣ ⎦ J . (4.9) The Jacobian matrix in Equation (4.9) relates the drive velocities, ( )tq , and accelera- tions, ( )tq , to the tool pose velocity and acceleration as: N N N N N ( ) , ( ) ( , ) x x y y z z a a aa a c c cc c x x xP P y y yP P z z zP P θ θ θθ θ θ θ θθ θ ⎡ ⎤ ⎡ ⎤⎡ ⎤ ⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥= = +⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦ ⎣ ⎦⎣ ⎦ ⎣ ⎦ x xQ Q Q J q J q J q q                          . (4.10) The simplified linear dynamics of a typical feed drive system is presented in Figure 4.2. The corresponding rigid body model can then be expressed as: Figure 4.2: Linear rigid body feed drive dynamics. Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 61 N N 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 xx x yy y zz z aa a c c c um cx x um cy y um cz z um ca a m cc c u ⎡ ⎤⎡ ⎤ ⎡ ⎤⎡ ⎤ ⎡ ⎤ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥+ = ⎢⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦⎣ ⎦ ⎣ ⎦ ⎣ ⎦ q qM C             N N x y z a c d d d d d ⎡ ⎤⎢ ⎥⎢ ⎥⎥ ⎢ ⎥−⎥ ⎢ ⎥⎥ ⎢ ⎥⎥ ⎢ ⎥⎣ ⎦ u d (4.11) where gitiaiii rKKJm = is the equivalent inertia, gitiaiii rKKBc = is viscous damping, )(tui is control voltage produced by the digital control law of the drive, )(tdi is equiva- lent disturbance, and the suffix , , , ,i x y z a c= denotes the drives, respectively. The amplifier and motor torque constants are, aK , tK , respectively, and gr is the gear reduc- tion ratio of the drives. The drive accelerations q can be represented in compact form by ( ) [ ( ) ( ) ( )]t t t t−= − −1q M u d Cq  . (4.12) By substituting the kinematics relationships from Equation (4.10) into Equation (4.12), the drive dynamics can be mapped to the tool pose ( )[ , , , , ]Tx y zP P P a c=x as: ( ) ( ) ( ) ( ) [ ( ) ( ) ( )]t t t t t t t−= + − −1x J q J M u d Cq   . (4.13) The actual position, velocity and acceleration of the pose can be written as a sum of the reference commands and corresponding tracking errors, , , , ( ) [ , , , , ] ( ) [ , , , , ] T x ref y ref z ref ref ref T Px Py Pz a c t P P P a c t e e e e e ⎫= ⎪ → = −⎬= ⎪⎭ ref ref x e x x e (4.14) Substituting Equation (4.14) into Equation (4.13), the tracking error dynamics of the tool is expressed by ( ) ( ) ( ) ( ) ( ) [ ( ) ( ) ( )]t t t t t t t t−= − − − −1refe x J q J M u d Cq   . (4.15) Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 62 Note that in Equation (4.15), the vector )(trefx is trajectory dependent. The Jacobian, )(tJ , and its derivative, )(tJ , introduce time varying non-linear coupling, since they depend on the drive position as well as the velocities. 4.3 Tool Tip Contour Error Estimation The reference trajectory of the tool tip, , , ,[ , , ] T x ref y ref z refP P P=refP , and corresponding tracking errors, [ , , ]TPx Py Pze e e=e , as well as the contour error vector, εG are illustrated in Figure 4.3a. The mutually orthogonal unit tangent, )(tt G , and normal, )(tnG , vectors along the path are computed using the reference velocity and acceleration of the tool tip: , , , , , ,, , , , , T T x ref y ref z ref x ref y ref z ref ref ref ref ref ref ref P P P P P P⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥= =⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦ t n P P P P P P      G G       (4.16) t n b F F e e h ε et ε F en+eben+eb Workpiece Coordinate System xy z xy z Desired Contour P xr ef , P yr ef , P zr ef P x,P y ,P z Actual Tool Tip Location Actual Tool Tip Location a) b) Figure 4.3: Modeling of tool tip contour errors in 5-axis machining. Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 63 where d dt=ref refP P and 22d dt=ref refP P . By defining the bi-normal vector )()()( ttt ntb GGG ×= , a moving local coordinate frame, namely, a Frenet Frame [70] on the tool-path at the reference tool center point ))(( trefP is defined as illustrated in Figure 4.3a: ( ) x x x y y y z z z t n b t t n b t n b ⎡ ⎤⎢ ⎥= ⎢ ⎥⎢ ⎥⎣ ⎦ F , (4.17) and the Frenet Frame in Equation (4.17) is augmented to the following transformation: 0 0 0 0 ( ) 0 0 0 0 0 1 0 0 0 0 0 1 x x x y y y z z z t n b t n b t t n b ⎡ ⎤⎢ ⎥⎢ ⎥⎢ ⎥= ⎢ ⎥⎢ ⎥⎢ ⎥⎣ ⎦ F (4.18) to include the rotary drive contributions. Hence, the following coordinate transformation holds: ( ) ( ) ( ) ( ) ( ) ( )t t t t t t= ↔ = TF Fe F e e F e , (4.19) where 1T FF −= and [ , , , , ] , ( ) [ , , , , ]T TPx Py Pz a c t n b a ce e e e e t e e e e e= =Fe e contain the tracking errors of the tool tip in the Workpiece Coordinate Frame and in the Frenet Frame, respectively. As shown in Figure 4.3a, the vector sums of the errors in normal and bi- normal directions represent an approximation to the contour error vector as: n be e≅ +εG G G . (4.20) Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 64 If the tool path does not contain sharp curvatures and the feeds are low, the contouring error vector expressed in Equation (4.20) is sufficiently accurate. However, at high speed contouring, the limited control bandwidth of the drives creates a time delay between the reference and actual positions of the tool. The tracking errors of the rotary drives may offset the workpiece and can even cause the actual tool position to lead the reference trajectory during 5-axis machining. The delay or lead time between the reference and actual tool tip positions is estimated from the error in tangent direction ( )te as: 2 2 2 , , , t d x ref y ref z ref et P P P = + +   , (4.21) and used to shift the Frenet Frame in time as follows: ( )dt t= −F F . (4.22) Since the reference trajectory is generated in advance, one can simply use the Shifted Frenet Frame ( )F and utilize the normal and bi-normal error components to improve the contour error estimation. Although the direction of the contour error vector is estimated more accurately with this method, such a transformation does not evaluate the amplitude of the tool tip contour error ( )εG accurately, leading to over or under-cuts during actual contouring control. As illustrated in Figure 4.3b, the vector between the current and the delayed reference positions is expressed as: 5 1 ( ) ( ) 0 , 0 dt t t × − −⎡ ⎤⎢ ⎥= ∈ℜ⎢ ⎥⎢ ⎥⎣ ⎦ ref refP P h h (4.23) Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 65 and its representation in F~ is computed by: = TFh F h . (4.24) In order to correct the estimation, only the normal and bi-normal components are used, and the tangent component is excluded by multiplying with (0, 1, 1,1,1)diag= − −W , =F Fh Wh (4.25) The tracking errors on the Compensated Frenet Frame, F~ become: = +TF Fe F e h  . (4.26) As shown in Figure 4.3b, the sum of vector components in the normal and bi-normal directions ( ˜ en + ˜ eb) evaluated from Equation (4.26) more accurately approximate the actual amplitude of contour error along arbitrary tool-paths as: 2 2n be e≅ +ε . (4.27) 4.4 Tool Orientation Contour Error Estimation As shown in Figure 4.4, there is a difference in definition of the regular tool orienta- tion tracking error vector oe G , and the orientation contour error vector ,actoεG , measured as the normal angular deviation from the orientation in the spherical coordinate system. The orientation tracking error vector is measured from the current tool orientation, O , to the current reference orientation. As shown by Equation (4.3), by simply minimizing the tracking errors of the rotary drives, orientation tracking errors can be eliminated. How- ever, doing so does not necessarily rotate the tool closer to the desired orientation Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 66 trajectory, refO . The tool can be rotated to the desired orientation trajectory by eliminat- ing the approximated contour error vector, ,act≅o oε εG G , as shown in Figure 4.4. Similar to the Cartesian contour error approximation, the orientation contour error vector is predicted by transforming the tracking errors to obtain the vector that is normal to the reference trajectory as: , , , o i o o j o o k ε ε ε ⎡ ⎤ ⋅⎢ ⎥= = −⎢ ⎥⎢ ⎥⎣ ⎦ o ω e ωε e ω ω G G GGG G G (4.28) where [ , , ]Ti j kω ω ω=ωG is the angular velocity of the tool. Normalizing the velocity of the tool axis, =ω ω ωG , Equation (4.28) can be re-written as: = − ⋅o o oε e ω(e ω) . (4.29) Figure 4.4: Modeling of tool orientation contour errors in 5-axis machining. Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 67 From the kinematics of the 5-axis machine tool, it is known that the orientation Jacobian ( )oJ in Equation (4.3) relates the rotary axis velocities )]v,v[( Tca=v to the orientation velocity of the tool through: k= oω J v , (4.30) where =v v v is the unit vector velocity in joint coordinates, and k = ω v is a scaling constant. In addition, assuming that the orientation errors are compensated with small drive rotations within the control sampling interval ( )1[ sec]T m= , we can use the following approximation: ≈o o Re J e (4.31) where [ , ]Ta ce e=Re is the array containing the individual rotary drive errors. Substitut- ing Equations (4.30) and (4.31) into Equation (4.29) yields: 2 T Tk= −o o R o o o Rε J e J vv J J e . (4.32) Using Equation (4.30) and properties of the unit vectors ( )ω, v , it can be shown that 2 21T T T T T T Tk k= = = → =o o o ov v ω ω v J J v v J J v , (4.33) which allows us to simplify Equation (4.32) as: T= − R o o R R ε ε J [e v(v e )] . (4.34) As a result, the relationship between the tool orientation contour error ( )oε in WCP and the corresponding rotary drive errors ( )[ , ]a cε ε=Rε in joint coordinates is obtained from Equation (4.34) as: Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 68 a T c ε ε ⎡ ⎤= = −⎢ ⎥⎣ ⎦R R R ε e v(v e ) . (4.35) 4.5 Implementation and Experimental Results The proposed contouring error prediction models are experimentally validated on a 5-axis machine tool, as modeled in Chapter 3. The in-house developed research CNC has an open architecture [72], which allows rapid implementation of trajectory generation and control laws. The Cartesian drives have 2.9 [micron] encoder resolutions, and the rotary table driven by the A and C drives are equipped with 8.7 [micro radian] resolution encoders. The dynamic parameters of the drives identified experimentally [47] are presented in Table 4.1. All five drives are individually controlled by robust sliding mode Table 4.1: Rigid body drive parameters. X-axis Y-axis Z-axis A-axis C-axis a t g Jm K K r = 0.00162 2Vsec mm⎡ ⎤⎣ ⎦ 0.00174 2Vsec mm⎡ ⎤⎣ ⎦ 0.00296 2Vsec mm⎡ ⎤⎣ ⎦ 0.00682 2Vsec rad⎡ ⎤⎣ ⎦ 0.00054 2Vsec rad⎡ ⎤⎣ ⎦ a t g Bc K K r = 0.00681 [ ]Vsec mm 0.00863 [ ]Vsec mm 0.01518 [ ]Vsec rad 0.01964 [ ]Vsec rad 0.00368 [ ]Vsec rad Table 4.2: Five-axis SMC control parameters (SISO scheme). SMC Parameters X-axis Y-axis Z-axis A-axis C-axis sec]/[radnω 200 200 180 150 200 sec]//[ mmVKs ([V/rad/sec] for rotary drives) 7 5 5 0.5 1 Disturbance Adaptation Gain [ / ]V mmΓ ( [ / ]V radΓ for rotary drives) 25 30 15 3 3 Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 69 controllers (SMCs) in order to obtain higher tracking performance. The design procedure is presented by Altintas and Erkorkmaz [35] where the sliding surface is based purely on the tracking errors of the individual axes. Thus, the SMC of each axis drive is tuned separately, and the bandwidth of each Cartesian axis is closely matched in order to obtain better contouring performance (see Table 4.2). The proposed contour estimation, as well as the robust axis control algorithms were implemented on the Labview RT© system comprised running on personal computers for real time computation of the control law, and a FPGA card for data acquisition and decoding of the encoder signals. The sampling rate of the control loop is fixed to 1kHZ. -30 -20 -10 0 -10 0 10 20 0 20 40 60 80 100 120 140 160 180 z [m m ] z [m m ] x [mm] i jy [mm] -1 1 1 -1 1 k a) b) Figure 4.5: Test Tool-path: a) reference tool tip ( ), , ,, , Tref x ref y ref z refP P P P⎡ ⎤= ⎣ ⎦ position, b) reference tool axis orientation ( ), , ,, , Tref i ref j ref k refO O O O⎡ ⎤= ⎣ ⎦ Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 70 The effectiveness of the proposed 5-axis contour error estimation algorithms is demonstrated on a tool-path that contains a circular contour for the tool tip motion and smoothly varying orientation for the tool axis. The reference 5-axis test tool-path is analytically given by 1 15cos( ) 5sin( ) [ ] 15sin( ) 5cos( ) [ ] 15 [ ] 0...4 , cos (sin( ) / 5) [ ] tan(cos( ) / 5) [ ] xref yref zref ref ref P l l mm P l l mm P l mm l a l rad c l rad π − ⎫= + ⎪⎪= − ⎪⎪= =⎬⎪= − ⎪⎪= − ⎪⎭ (4.36) and shown in Figure 4.5. The reference path of the tool orientation continuously changes as shown in Figure 4.5b. The path of the tool tip follows a circular motion in the x-y plane, and at the same time the cutter moves in the z direction, where the total distance of travel along the workpiece is about 273.87 [mm]. The reference tool-path is interpo- lated by generating a smooth trajectory with a maximum cruise feed of 50 [mm/sec] and cubic acceleration/decelerations not exceeding 1500 [mm/sec2]. Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 71 In order to investigate the accuracy of the proposed 5-axis tool tip contour error computation methods, the basic Frenet Frame (Equation (4.20)), Shifted Frenet Frame (Equation (4.22)) and the Compensated Frenet Frame (Equation (4.26)) are compared against true tool tip contour errors as shown in Figure 4.6. The true contour errors are calculated through an iterative search of the shortest distance between the tool tip position and actual path, as presented by Altintas and Erkorkmaz [39], measured in the workpiece coordinate system. The proposed Compensated Frenet Frame approximation provides the most accurate estimation of contour errors as shown in Figure 4.6a. Severe peaks in the contour error trend can be observed at times ~1, 2.2, 3.5 and 5 [sec] during 0 200 400 a) True Tool Tip Contour Error(solid) & Proposed Estimation(dashed) C on to ur E rr or [m ic ro n] D is cr ep en cy [m ic ro n] 0 0.2 0.4 0.6 b) 0 1 4 5 1 2 4 5 1 2 2 4 5 -5 0 5 time [sec] time [sec] time [sec] [m ili se co nd ] c) Compensated Frenet Frame Shifted Frame Basic Frame Figure 4.6: Comparison of tool tip contouring error estimations (experimental): a) com- parison of true and estimated tool tip contour errors, b) discrepancy of proposed estimation from true contour errors, c) delay time Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 72 air cutting. In addition, relatively smaller jumps occur in between at ~1.5, 3, 3.5 and 4.2 [sec]. As the rotational axes go through a direction reversal during 5-axis operation, the joint velocity changes rapidly with high acceleration, and, at the same time, the coulomb friction at the joint acts as a sudden step disturbance to the drive. Each axis controller deals with the disturbance individually and tries to track the high acceleration reference command. The disturbance is rejected in finite time, and severe tracking errors at the axis level occur, jeopardizing the coordination between the drives. This deteriorates the contouring performance severely. Thus, the contour error peaks in Figure 4.6a coincide with the locations when the A rotary axis goes through velocity reversal, and the smaller error peaks are related to the tracking performance of the C axis (see Figure 4.7a and b). This indicates that during 5-axis machining not only the Cartesian but also the rotary drives contribute strongly to the contouring performance. A closer view of the estimation errors of the basic, delayed and compensated Frenet Frames is shown in Figure 4.6b. All of the proposed methods estimate the contour errors within an acceptable average estimation discrepancy (0.01~0.05[micron]). The most severe cases occur when the contour error spikes above 300 [micron] where the actual tool tip has its maximum lag from the reference point. This lag can be observed from the delay time computed from Equation (4.22) and plotted in Figure 4.6c. The basic and Compensated Frenet Frame approximations show good performance even at these locations with an estimation error as high as 0.1 and 0.05[micron], respectively. As a result, both the basic and the Com- pensated Frenet Frame methods prove to be sufficiently accurate for contouring control even when the tool tip deviates highly from the desired path. On the other hand, the delayed Frenet Frame shows the worst performance confirming that the offset h in Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 73 Equation (4.23), must be compensated to estimate the amplitude of the contour error accurately. Similar to the Cartesian case, the actual orientation contour error vector can be ob- tained by iteratively finding the closest normal point to the actual reference orientation tool-path in spherical coordinates. The actual orientation contour error vector ( )o,actε is computed and transformed into the corresponding a and c rotary drive tracking errors ,( )a cε using the orientation Jacobian (Equation (4.3)): ,a c ≈ Oε J ε , (4.37) 0 1 4 5 −3 −2 −1 0 1 2 3 a) b) c) d) time [sec] Er ro r [ m ill ira di an s] 0 1 4 5time [sec] 0 1 4 5time [sec] 0 1 4 5time [sec] −2 −1 0 1 2 0 0.001 0.002 -0.002 -0.001Es tim at io n Er ro rs [m ill ira di an s] 0 0.001 0.002 -0.002 -0.001 Tracking Errors Tracking Errors True&Approximated Contour Error Components True&Approximated Contour Error Components Figure 4.7: Comparison of tool orientation contour error estimations: a) A axis errors, b) C axis errors, c) discrepancy in the proposed contour error estimation in A axis, d) discrepancy in the proposed contour error estimation in C axis. Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 74 and compared against the proposed analytical estimation given in Equation (4.35). Figure 4.7 represents the experimentally recorded rotary drive tracking errors, as well as the true and estimated orientation contour errors expressed in drive coordinates. As shown in Figure 4.7a and b, the transformation of the orientation contour errors com- puted from Equation (4.31) and (4.35) match, and they both show distinct error trends in comparison to the tracking errors of the rotary drives. Rotary drive tracking errors show peaks due to the coulomb friction at path locations when A and C axes go through velocity reversals. The orientation contour error component mainly follows the tracking errors, indicating that the motion should be synchronized. However, deviations between the trends occur especially when synchronization is lost, as one rotary drive is exposed to a velocity reversal yielding high tracking errors and the other one is not. In Figure 4.7a, one can see that the A rotary drive tracking and orientation error trends differ from each other at 0.5, 1.5, 3, 4 [sec] on the time scale, corresponding to the locations where the A rotary drive has its highest tracking errors. Similar observations can be made for the C axis at 1, 2.5, 3.5 and 5 [sec], indicating that elimination of orientation contour errors is crucial in establishing synchronization of the drives towards improving the part accuracy. Figure 4.7c-d show the deviation of the proposed analytical orientation con- tour error estimation from the measured values. The maximum deviation is less than 6 [microradian] for both A and C drives. Considering that the rotary encoders have a resolution of 8.7 [microradian], the proposed estimation is valid and suitable for design- ing orientation contour controller that is presented in Chapter 5 of this thesis. Chapter 4: Modeling of Contouring Error for 5-Axis Machine Tools 75 4.6 Conclusion The main concerns during 5-axis contouring are the tool tip and the tool orientation contour errors. The kinematic model of 5-axis machine tools shows that all five physical drives influence tool tip positioning accuracy, while the tool axis orientation errors are caused only by the rotary drives. In this chapter, the kinematics and rigid body dynam- ics of five axes, namely, three translational and two rotary drives are combined to predict tool tip and tool orientation errors that produce contouring errors on the machined part surface. It is shown that even if the tracking accuracies of all five drives are individually greatly improved by high bandwidth, robust controllers, the overall contour error along 5-axis tool paths can still be significant due to the path curvature and the kinematics of the machine. Instead of controlling each drive individually, the experimentally vali- dated, real-time, contouring error estimation model allows the design of a global, sliding mode controller which simultaneously controls all five drives with the objective of achieving high path tracking accuracy. The contributions of this chapter have been published in [26], and a contour controller that improves the path tracking accuracy is presented in Chapter 5 of this thesis. 76 5. Precision Contouring Controller Design for 5-Axis Machine Tools 5.1 Introduction Two types of contour errors were introduced in Chapter 4. The first type of contour error is the normal deviation of the tool tip from the desired tool-path, called the tool tip contour error, arises as a nonlinear function of the tracking errors of all five drives. Orientation contour error, on the other hand, has been defined as normal deviation from the reference tool orientation trajectory, and it can be controlled only by the rotary drives of a 5-axis machine tool. Design of a multi-input-multi-output sliding mode controller for minimizing both contour errors is introduced in this chapter of the thesis. The proposed 5-axis contour controller is designed in two steps. In Section 5.2, the kinematics of the machine tool and the tool tip contour error models introduced in Chapter 4 are utilized to design a robust sliding mode controller that minimizes the tool tip contour errors. The proposed tool tip contour error controller is based on the rigid body dynamic model. It estimates the normal deviation of the tool tip from the desired contour and generates coupled control to the axis drives (X,Y,Z, A and C). The design of an orientation contour error controller is presented in Section 5.3. The orientation sliding mode controller, on the other hand, utilizes the orientation contour error estimation and couples only the two rotary drives (A and C) for accurate orientation contouring. Both contour controllers are implemented together, and the complete 5-axis contour controller scheme is constructed to control the tool tip and orientation contour errors in a coupled fashion. The effectiveness of this control strategy is experimentally demonstrated on the 5-axis CNC machine tool in Section 5.4, and the contributions [43] are summarized in Section 5.5. Chapter 5: Precision Contouring Controller Design for 5-Axis Machine Tools 77 5.2 Design of the Sliding Mode Controller for Tool Tip Contour Errors The simplified rigid body dynamics of a typical feed drive system on a 5-axis ma- chine tool were modeled in Chapter 4 as: ( ) [ ( ) ( ) ( )]t t t t−= − −1q M u d Cq  (5.1) where ( ) [ ( ), ( ), ( ), ( ), ( )]Tt x t y t z t a t c t=q contains the drive positions. q and q contain the drive velocities and accelerations, respectively. u is the control input to the amplifiers, and d is the external disturbance reflected at the amplifier’s input. ∈ℜ5×5M and ∈ℜ5×5C are diagonal matrices that contain the drives’ equivalent inertia and viscous damping values. The jacobian matrix, ( )t ∈ℜ5×5J relates drive velocities, q , to the tool pose velocity ( ( ) [ ( ), ( ), ( ), ( ), ( )] )Tx y zt P t P t P t a t c t=x      and obtained from the kinematics of the 5-axis machine. The drive dynamics are mapped to the tool pose using the jacobian, and tracking error dynamics of the tool in the workpiece coordinate system has been expressed as: ( ) ( ) ( ) ( ) ( ) [ ( ) ( ) ( )]t t t t t t t t−= − − − −1refe x J q J M u d Cq   . (5.2) Based on the contour error and drive dynamics models presented in Chapter 4, a sliding mode controller that minimizes the tool tip contour errors is designed herein. Chapter 5: Precision Contouring Controller Design for 5-Axis Machine Tools 78 The actual tool tip contour error vector ( )ε has been modeled by reflecting the tool tracking tip errors ( )e on the delay compensated Frenet Frame ( )F , ⎫= + ⎪⎪= + + ⎬⎪= + + + ⎪⎭ T F F T T F F T T F F e F e h e F e F e h e F e 2Fe F e h                  . (5.3) The contour errors left on the part surface are estimated from the normal and bi-normal components ( )n be eε = +G   of the errors (see Figure 5.1), which are minimized by the sliding mode controller. The contour error dynamics are derived by substituting Equa- tion (5.2) into Equation (5.3): F e h ε F en+eb Workpiece Coordinate System xy z Actual Tool Tip Location Figure 5.1: Tool tip contour error estimation. Chapter 5: Precision Contouring Controller Design for 5-Axis Machine Tools 79 ( ) t n b a c e e e e e − ⎡ ⎤⎢ ⎥⎢ ⎥⎢ ⎥ ⎡ ⎤= = − − + − + + +⎣ ⎦⎢ ⎥⎢ ⎥⎢ ⎥⎣ ⎦ T 1 T T F ref Fe F x Jq JM u d Cq 2F e F e h                . (5.4) Note that the dynamics of the first three components, namely the tangent, normal and bi- normal ( ), ,t n be e e   errors are nonlinear and time varying due to the Jacobian and Frenet Frame transformations. The last two components are the regular tracking error dynamics of the rotary drives. A second order sliding surface, 5 1×∈ℜS , is designed to contain the proportional, integral and derivative of the time varying errors according to: 0 t n t b a c S S S d S S τ ⎡ ⎤⎢ ⎥⎢ ⎥⎢ ⎥= = + +⎢ ⎥⎢ ⎥⎢ ⎥⎣ ⎦ ∫P F I F D FS C e C e C e   , (5.5) where 55×ℜ∈PC , 55×ℜ∈IC and 55×ℜ∈DC are diagonal design matrices that represent the desired achievable dynamics of the errors on the sliding surface. The control input ( )5 1×∈ℜu must be manipulated in such a way that the errors and the time derivative converge asymptotically to the stable sliding surface, so that they eventually slide to origin →S 0 as 0eF →~ and 0eF →~ . The error dynamics on the sliding surface can be represented in Laplace domain as: 2 ( ) 0s s s⎡ ⎤= + + =⎣ ⎦D P I FS C C C e  (5.6) Chapter 5: Precision Contouring Controller Design for 5-Axis Machine Tools 80 where s is the Laplace operator, and stable 2nd order error dynamics can be designed by choosing constant design matrices: ( ) ( ) 5 5 , , , , , 2 2 2 2 2 , , , , , 2 ,2 ,2 ,2 ,2 , , , , t n t n n n b n b a n a c n c n t n n n b n a n c diag diag ζ ω ζ ω ζ ω ζ ω ζ ω ω ω ω ω ω × ⎫= ⎪⎪= ⎬⎪= ⎪⎭ D P I C I C C (5.7) where iζ and in,ω , cabnti ,,,,= are the individual damping and natural frequencies for each tracking error component. The variation in the external disturbance caused by the cutting process and friction is considered to be constant and remains between the known upper, 5 1+ ×∈ℜd , and lower, 5 1− ×∈ℜd , limits during short time intervals. The external disturbances force contouring errors to deviate from the sliding surface, but they can be tracked using the following observer which integrates the sliding surface: ( )ˆ ˆ ˆ 1k T= → = − +T Td ΓκJ FS d(k) d ΓκJ FS   . (5.8) T is the control sampling period, ( ), ,...,t n cdiag= Γ Γ ΓΓ is the observer gain matrix, TJ is the transpose of the jacobian, and ( ), ,...t n cdiag κ κ κ=κ is a matrix with positive entries used to impose limits on the integral control action against the disturbances, so that the observations are within the given bounds, ˆ− +≤ ≤d d d , ˆ0 if and ˆ0 if and 1 otherwise − + ⎧ ⎫≤ ≤⎪ ⎪⎪ ⎪= ≥ ≥⎨ ⎬⎪ ⎪⎪ ⎪⎩ ⎭ d d S 0 κ d d S 0 . (5.9) Chapter 5: Precision Contouring Controller Design for 5-Axis Machine Tools 81 The bounded disturbance observer avoids chattering of the controller on the sliding surface. In order to push the tracking errors onto the sliding surface, the following Lyapunov vector function is postulated as: 1 ˆ ˆ 2 t n b a c V V V V V − ⎡ ⎤⎢ ⎥⎢ ⎥ ⎡ ⎤⎢ ⎥= = + − −⎣ ⎦⎢ ⎥⎢ ⎥⎢ ⎥⎣ ⎦ T T 1V S MS (d d) Γ (d d) . (5.10) This function penalizes the deviation from the surface and the square of the disturbance estimation error. The derivative of the Lyapunov function, ( )ˆ ˆ−= − − TT 1V S MS d d Γ d (5.11) must be negative for asymptotic stability. This pushes the states towards the sliding surface where they follow the desired 2nd order dynamics. Equations (5.6) and (5.8) are substituted into Equation (5.11) to yield: ( )ˆ⎡ ⎤= + + − −⎣ ⎦ TT TI F P F DV S M C e C e C e d d κJ FS     (5.12) and Equation (5.4) is used to expand the derivative of the Lyapunov function: ( )ˆ − −⎡ ⎤= + + − − + + + +⎣ ⎦ + − − T T T T 1 T 1 T T I F P F ref T T T T V S M C e C e F x F Jq F JM u F JM Cq 2F e F e h S F Jd S F J d d            . (5.13) Noting that ( )ˆ− −T T T TS κF Jd S F J d d in Equation (5.13) can be regrouped as ( )( )ˆ− −T TS F J d d 1 κ , and due to the disturbance limits imposed in Equation (5.9) Chapter 5: Precision Contouring Controller Design for 5-Axis Machine Tools 82 ( )( )ˆ− − ≤T TS F J d d 1 κ 0 holds. As a result, the following criterion is sufficient to ensure the asymptotic stability ( )( )t , ,min 0nω > , and max min 0ζ ζ> > . It should be noted that upper bound on maxζ is used to avoid the 3rd pole ( )23 1n np jζω ω ζ= − − − from approaching zero, and the allowable pole locations are illustrated in Figure 6.4. Hence, the simple Least Squares problem in Equation (6.10) is constrained by Equation (6.14) as: ( ) ( )1min subject to: 0 2 ≥T θ Y -Φθ Y -Φθ h . (6.15) The augmented problem in Equation (6.15) is a Linear Programming problem with non- linear constraints has been solved using the Particle Swarm Optimization method as described in the following section. Chapter 6: Closed Loop Identification of Feed Drive System 112 6.4.1 The Principle of PSO Particle Swarm Optimization (PSO) algorithm [62] is a population based intelligent search where individuals called as ‘particles’, ,1 ,2 ,[ , , ] T i i i i Dσ σ σ=σ … containing the optimization variables (1... )D that are essentially candidate solutions for a given optimi- zation problem. Each particle’s ‘position’, i.e. the optimization variables assigned to that particle, is adjusted dynamically by the amount of ,1 ,2 ,[ , , ] T i i i i Dv v v=v … , denoted as the particle’s velocity. The PSO algorithm keeps track of each particle’s history, and stores its best value, ,i pbestσ that has minimized the objective function. Another best value that is tracked by the PSO minimizer is the overall best value of the entire particle batch obtained so far; ,min( ), 1...gbest i pbest i M= =σ σ , where M is the total number of particles. Hence, the particles’ position and velocity are adjusted each iteration l with respect its own best location and the best location of the group as: ( ) ( )1 1 , 2 1 1... l l l l l l i i i pbest i gbest i l l l i i i w i M η η+ + = ⋅ + − + − = + = v v σ σ σ σ σ σ v . (6.16) Im Re Allowable Root Locations Im Re −ωn,min −ωn,min −ζn,min −ζn,max −pmin p1 p2 p3 −ωn,min(ζmax+ ζmax-1) a) b) Figure 6.4: Allowable locations for; a) real, b) complex poles. Chapter 6: Closed Loop Identification of Feed Drive System 113 ( )1 2, 0,2η η ∈ are the learning gains that represent the weighting of the stochastic accele- ration terms, which pull each particle towards ,i pbestσ and gbestσ . Low learning values allow particles to search far from the target, and high values cause quick convergence with abrupt movement toward and past the optimum. w is the inertia weight, and it controls the convergence behavior of the PSO [74], [75]. 6.4.2 PSO Solution In setting up the PSO to solve the optimization problem given in Equation (6.15), the solution search is conducted directly in terms of the closed loop pole frequency ( nω , p ) and damping ( )ζ values. As shown in Figure 6.5, the search region is bounded by the constraints expressed in Equation (6.14); thus, adopting this strategy simplifies the incorporation of stability constraints. Each particle contains the optimization va- riables, namely the pole locations ( nω , p ) and the damping value ( )ζ : Figure 6.5: Search region for the Particle Swarm Optimization algorithm. Chapter 6: Closed Loop Identification of Feed Drive System 114 , , 1... T i i n i ip i Mω ζ⎡ ⎤= =⎣ ⎦σ . (6.17) As the particles sweep the bounded search region, the objective function, ( ) ( )1 2 f = TY -Φθ Y -Φθ (6.18) is evaluated. Therefore, the parameter vector, θ is separated into two parts: [ ] 1 2 3 1 2 1 2 1 2[ ] TT oα α α β β β δ= = θ θ θ θ θ  , (6.19) and the particles (i.e. optimization variables) are related to the first part of the drive parameters [ ]1 3 1 2 Tα α α=θ by using Equations (6.8) and (6.13) as: 2 3 1 3 1 22 2 2 2 2 2 2 1 1, , 2 2 2 n n n n n n n n a p a p p a p a p a p ω ωα α αω ω ω ω ω ω += = = = = =+ + + . (6.20) In order to relate the particles to the remaining drive parameters, [ ]2 0 1 2 Tβ β β δ=θ , the Least Squares sub-problem [65] is constructed from Equation (6.10). Upon the existence of ( ) 12 2 −TΦ Φ , 2θ is obtained uniquely from: [ ] ( ) ( )111 2 2 2 2 2 1 1 2 −⎡ ⎤ = ⇒ = −⎢ ⎥⎣ ⎦ T TθΦ Φ Y θ Φ Φ Φ Y Φ θ θ . (6.21) By using only the pole frequency ( nω , p ) and damping ( )ζ values as particle posi- tions ,( [ ] ) T i i n i ip ω ζ=σ , the objective function ( ) ( )( )1 2f = TY -Φθ Y -Φθ is minimized by implementing the PSO algorithm as follows: Chapter 6: Closed Loop Identification of Feed Drive System 115 Step 1. Specify the search region boundaries of the particles as min min ,min min T npσ ω ζ⎡ ⎤= ⎣ ⎦ and max ,max max Tnpσ ω ζ⎡ ⎤= ⎣ ⎦ according to the constraints, and initialize each particle’s position and velocity randomly as: 0 max 1 max 0 min 2 max min 2 ( ) i i v v rv rσ σ σ σ ⎫= − + ⎪⎬= + − ⎪⎭ (6.22) where r1 and r2 are random numbers between 0 and 1 where max max min( ) (10% ~ 20%)v σ σ= − × Step 2. Evaluate the objective function lif for each particle using Equations (6.18), (6.20) and (6.21). Step 3. Find , l i pbestf and compare it against lif for each particle. If , l l i i pbestf f≤ then set , l l i best if f= , and , li pbest i=σ σ . Step 4. Search for the minimum value of lif among all the particles; denote it as lgbestf , and set l gbest i=σ σ . Step 5. Check each particle individually, and if any value of nω , p and ζ is outside the search region, it is replaced within the closes bound as: min min max max , ,min , ,min , ,max , ,max min min max max if , if if , if if , if l l l l i i i i l l l l n i n n i n n i n n i n l l l l i i i i p p p p p p p p ω ω ω ω ω ω ω ω ζ ζ ζ ζ ζ ζ ζ ζ < ⇒ = > ⇒ = < ⇒ = > ⇒ = < ⇒ = > ⇒ = . (6.23) Step 6. Return to step 2 and iterate for the next particle position and velocity ( 1li +σ and 1li +v ) using Equation (6.16). Chapter 6: Closed Loop Identification of Feed Drive System 116 Step 7. When the number of iterations reaches the maximum, the search algorithm is terminated, and the latest gbestσ containing the optimal pole locations and damping values for the model is accepted. The model parameters ( )θ are computed from Equa- tions (6.20) and (6.21) and de-normalized by Equation (6.8) to be used in Equation (6.3). 6.4.3 PSO Search Parameter Tuning Although PSO is a stochastic algorithm, some guidance is required for selecting proper parameters to ensure convergence. The learning gains 1 2,η η and the inertia weight w control the performance and convergence of the PSO algorithm [76]. Assum- ing that each particle contains only one element, D=1 and during steady state the ,i gbestσ and ,i pbestσ values are constant, Equation (6.16) can be re-written for each particle as: N N ( ) 1 1 1 , 2 1 21 1 , and 1 1 ll l l l li i i pbest gbestl l i i u wv v w η ησ η σ η η ηησ σ + + + ⎡ ⎤ ⎡ ⎤−⎡ ⎤ ⎡ ⎤= ⋅ + + = +⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥−⎣ ⎦ ⎣ ⎦⎣ ⎦ ⎣ ⎦ BH xx    . (6.24) Equation (6.24) models the evaluation of the particles as a discrete linear dynamic system in state space form. [ ]l l l Ti iv σ=x is the state vector, [1 1]T=B is the input matrix and 1 , 2 l l i pbest gbestu ησ η σ= + is the constant input. The stability of the discrete PSO model can be ensured by setting the eigenvalues of the state matrix H within the unit circle [28]. The characteristic equation is calculated as ( ) 2det ( 1)w wλ λ η λ− = + − − +I H , (6.25) and the roots are: Chapter 6: Closed Loop Identification of Feed Drive System 117 2 21,2 [( 1 ) ] 2, 0 , ( 1) 2 ( 1) [( 1 ) ] 2, 0 w w w w j ηλ η ηη ⎧ ⎫+ − ± ∆ ∆ ≥⎪ ⎪= ∆ = − − + +⎨ ⎬+ − ± −∆ ∆ <⎪ ⎪⎩ ⎭ . (6.26) The objective is just to ensure convergence of the PSO algorithm, not tuning of the convergence dynamics. Therefore, the following simple observations give a quick guideline in selecting the PSO parameters. According to Equation (6.26), when 0∆ < , 1,2λ are complex conjugate, their norm become 1 2 wλ λ= = . Therefore, setting (0,1)w∈ ensures the convergence of PSO. In parallel to this fact: 1. When 0η = , 1 21, wλ λ= = , the PSO is at the border of convergence or diver- gence. 2. When 1η = , 2 4 0w w∆ = − < , and for ( )0,1w∈ , the PSO converges. 3. When 2η = , 2 26 1 ( 3) 8w w w∆ = − + = − − and for ( )0.172,1w∈ , 0∆ < PSO converges. 4. When 3η = , 2 28 14 ( 4) 12w w w∆ = − + = − − and for ( )0.536,1w∈ , 0∆ < , the PSO converges. 5. When 4η = , 2 210 9 ( 5) 16 0w w w∆ = − + = − − > and 2 [( 3) ] 2 1wλ = − − ∆ < − , it results in divergence. As a result, the learning gains can be selected as ( )1 2 0,3η η η+ = ∈ , and the weight factor is determined in the interval of ( )0.6,1w∈ to ensure the convergence of PSO. More analysis can be performed in order to select the optimal parameters for quick convergence by imposing faster dynamics on the system. Chapter 6: Closed Loop Identification of Feed Drive System 118 6.5 Simulation and Experimental Results The effectiveness of the identification strategy has been verified in simulations and experimental case studies conducted on virtual and actual feed drives. Since the virtual axis dynamics are already known, simulations allow a comparison between the true and identified drive parameters using both LS and PSO methods. The identification G-code used in tests consists of 70 random linear movements commanded within the range of 50[mm]± with a maximum feedrate of 250[mm/sec] and 22500[mm/sec ] acceleration. Due to the short distance between the position com- mands, the desired velocity is not necessarily reached, which allows the drive response to be observed for a wider velocity range. The identification reference trajectory is shown in Figure 6.6. Trapezoidal velocity profiling is being used to ensure good persis- Figure 6.6: Reference identification trajectory send to the drives. Chapter 6: Closed Loop Identification of Feed Drive System 119 tence in the excitation signal, and the reference and actual position commands are collected at a sampling period of 1[msec]sT = . The PSO constrained identification technique was validated in virtual drive parameter estimations by comparing the results against basic LS algorithm that does not accommodate any pole constraints. In applying the PSO, the search space is bounded depending on the sampling frequency as ,min min 0.2[Hz]n pω = = , ,max max 500[Hz]n pω = = , ,min 0.2nζ = and ,max 2.0nζ = . The weights are set to 1 20.103, 2.897η η= = and 0.6w = . The number of particles is 200M = and the maximum number of iterations is set to 100. 6.5.1 Validation on a Virtual Drive The identification strategy was first evaluated on a P-PI cascade controlled drive system. The true, LS and PSO estimated closed loop parameters are summarized in Table 6.1. A comparison between the true and estimated tracking and disturbance frequency response functions (FRF) is presented in Figure 6.7. The LS identification gives the most accurate results in terms of minimizing the objective if there is not any stability constraint active. As can be noted from the table, there is some discrepancy between the true and LS estimated parameters. The slower pole at 10.8 [Hz] is estimated at 11.4 [Hz], and the complex conjugate poles at 30.6 [Hz] are estimated at 31.3 [Hz]. Table 6.1: Actual and estimated closed loop parameters for P-PI controlled virtual drive(simulation). Search Variables Actual LS Estimate PSO Estimate [Hz]nω 30.6 31.3 31.3 []ζ 0.45 0.39 0.39 p[Hz] 10.8 11.4 11.4 Chapter 6: Closed Loop Identification of Feed Drive System 120 Model Parame- ters Actual LS Estimate PSO Estimate a1 238.4335 226.9455 226.94681 a2 48421.5859 49915.08 49915.0545 a3 2.487140.30 2778133.78 2778130.868 b0 0 -0.05340 -0.05340 b1 117.5424 104.3654 104.366781 b2 37258.9796 37446.2493 37446.2295 dc 2.6815 585.30 585.31 Considering the FRFs shown in Figure 6.7a and b, the identified transfer functions are in good agreement up to the frequency range 30[Hz], which is within the tracking range of most CNC motion commands. The true friction value is obtained by averaging the positive and negative friction torques of the virtual drive. As observed, there is 218× mismatch between the estimated and true friction values. As can be noted from the phase match, disturbance transfer functions are in good agreement. There is 221× difference 10−1 100 101 a) M ag ni tu de [m m /m m ] 100 101 102 −150 −100 −50 0 50 Ph as e[ de g] 10−6 10−5 10−4 10−3 10−2 b) 100 101 102 −200 −150 −100 −50 0 50 100 Frequency[Hz] Frequency[Hz] LS and PSO Estimation Actual LS and PSO Estimation Actual Figure 6.7: True and estimated FRFs for P-PI controlled drive system (simula- tion): a) tracking transfer function, b) disturbance transfer function. Chapter 6: Closed Loop Identification of Feed Drive System 121 between the DC amplitudes, which actually compensates for the discrepancy between the true and estimated friction parameters. Overall, the LS estimation is successful in capturing the key dynamic characteristics required to predict the tracking performance of the actual drives. The PSO, on the other hand is slightly less accurate in minimizing the objective function, since it is not analytical. Figure 6.8 shows the convergence of the PSO estimated drive parameters. As noted, the parameters converge to the correct values with a small discrepancy well within 100 iterations. The actual and identified axis dynamics have also been compared in trajectory tracking tests to prove the practicality of the methods, and the results are presented in Figure 6.9. The true axis response agrees with the LS and PSO model closely proving that the PSO is accurate enough to be implemented for solving the optimization problem. 10 20 30 40 50 60 70 80 0 0.2 0.4 0.6 0.8 1 Objective Fucntion Value 50 0 0.3 200 400 600 800 ωn ζ [r ad /s ec ] 50 0.4 0.5 0.6 0.7 [] [] 50 40 50 60 70 80 p [r ad /s ec ] Iteration Iteration Iteration Iteration Figure 6.8: Parameter convergence of the PSO algorithm for P-PI controlled drive system. Chapter 6: Closed Loop Identification of Feed Drive System 122 The second virtual experiment is performed on a PID controlled system and the es- timated drive parameters are summarized in Table 6.2. As shown in the table, there is a close match between the actual and identified drive parameters that can also be verified from the FRFs presented in Figure 6.10. The tracking transfer functions are in good agreement up to 50 [Hz]. The PSO and LS estimates are almost identical, with slight discrepancy in nω and ζ . The tracking responses of the estimated drive models are compared in Figure 6.12, which again shows very close agreement between the true system response and both LS and PSO identified model responses. The discrepancy between the true and identified tracking errors is within 2 [micron]. Figure 6.9: Predicted and actual tracking performance of the P-PI controlled virtual feed drive system (simulation). Chapter 6: Closed Loop Identification of Feed Drive System 123 Table 6.2: Actual and estimated closed loop parameters for PID controlled virtual drive(simulation). Search Variables Actual LS Estimate PSO Estimate [Hz]nω 49.1 52.53 53.1 []ζ 0.707 0.62433 0.688 p[Hz] 7.4 7.63 7.62 Model Parameters Actual LS Estimate PSO Estimate a1 475.31142 460.073072 506.254774 a2 111794.743 128718.4765 133200.5673 a3 4263669.10 521019.0346 5327596.9951 b0 0 -0.13066504 -0.1758883 b1 286.861167 273.9933732 302.7291668 b2 76074.4037 85069.40516 88568.0002 dc 2.6815 626.397 646.906 150 100 10−1 100 101 a) M ag ni tu de [m m /m m ] 100 101 102 100 101 102 − − −50 0 50 Ph as e[ de g] 10-6 10-5 10-4 10-3 10-2 b) M ag ni tu de [m m /m m ] −200 −150 −100 −50 0 50 Ph as e[ de g] Frequency[Hz] Frequency[Hz] LS and PSO Estimation Actual LS and PSO Estimation Actual Figure 6.10: True and estimated FRFs for PID controlled drive system (simula- tion): a) tracking transfer function, b) disturbance transfer function. Chapter 6: Closed Loop Identification of Feed Drive System 124 6.5.2 Validation on a Ball Screw Drive After the identification algorithms have validated in virtual environment, they have been tested experimentally on an experimental ball screw drive setup to demonstrate practical effectiveness. Identification of a P-PI controlled drive system is first tested. A trapezoidal trajectory with 250[mm/sec] cruise velocity and 22500[mm/sec ] accelera- tion is used to validate the tracking estimation and tracking results are presented in Figure 6.12a. P-PI controlled true axis tracking errors are observed in the range of 400[micron] , and the identified model predicts the true tracking error trend with a maximum discrepancy of 20[micron]. The Sliding Mode Control (SMC) scheme has also been implemented in the ball screw drive setup. The SMC cancels out inertia and friction to widen up the position bandwidth, resulting in faster dynamics that are more difficult to excite with the reference trajectory for closed loop identification. Hence, it has been observed that the LS identification procedure has produced unstable poles locations in many tests due to the pole zero cancellations. The PSO algorithm on the other hand, has identified stable model parameters. As shown in Figure 6.12b, the Figure 6.11: Predicted and actual tracking performance of the PID controlled virtual drive system (simulation). Chapter 6: Closed Loop Identification of Feed Drive System 125 tracking errors for SMC system are in the range of 40 [micron], and some chattering is observed in the tracking response. The LS and PSO estimated models still capture the general drive response, and the average Coulomb friction estimates accurate enough allowing the identification method to be used in estimating systems with wider tracking bandwidths. Finally, the rapid identification technique has been tested on a virtual 5-axis ma- chine tool for simulating the contouring errors. P-PI cascade controllers with feedforward dynamics compensation are designed for each motion axis in order to achieve moderate tracking performance. The true and estimated drive parameters are compared in Table 6.3. Since the LS algorithm produces unstable poles with increased feedforward dynamic compensation, only the PSO is used to identify the drive parame- ters. As noted from the table, the X and Y poles are identified very closely. On the other −100 50 [m m ] Test Trajectory Test Trajectory a) b) −500 0 500 −40 40Tracking Errors (True,LS & PSO Estimations) Tracking Errors (True,LS & PSO Estimations) [m ic ro n] 0 0.5 1.5 2 2.5 0 0.5 1.5 2 20 Discrepency in Tracking Error Prediction [m ic ro n] Reference,Actual and Estimations True and Estimate True Estimated Time[sec] Time[sec] 10 Figure 6.12: Predicted and experimentally verified tracking performance of P-PI (a) and SMC (b) controlled ball screw drive: a) P-PI control results, b) SMC control results. Chapter 6: Closed Loop Identification of Feed Drive System 126 hand, the Z, A and C axis poles do not match. Because of the pole-zero cancellations, some of the poles could not be identified correctly. Table 6.3: Actual and estimated drive parameters of a P-PI controlled virtual 5-axis machine tool. Actual X Y Z A C Se ar ch V ar ia bl es [Hz]nω 55.54 73.5 73.6 31.9 27.3 []ζ 0.39 0.27 0.28 0.47 0.57 p[Hz] 19.3 22.12 22.2 9.29 8.45 Po le s a nd Z er os (F re q. [H z] /D am p. []) p1 55.7 0.39 73.5 0.27 73.6 0.27 31.9 0.47 27.3 0.57 p2 55.7 0.39 73.5 0.27 73.6 0.27 31.9 0.47 27.3 0.57 p3 19.3 1 22.1 1 22.2 1 9.29 1 8.45 1 z1 - - 156 0.72 72.4 -1 43.9 0.34 27.2 27.2 z2 46.1 1 156 0.72 72.4 -1 43.9 0.34 0.47 0.47 z3 25.6 1 24.4 1 22.7 1 9.88 1 8.61 1 dc 2.6815 1.835 183.57 0 0.2 PSO Estimated Se ar ch V ar ia bl es [Hz]nω 58 90.9 138.5 24.19 68.9 []ζ 0.27 0.2 2 1.53 2 p[Hz] 20.12 24.5 -138.5 24.19 68.9 Po le s a nd Z er os (F re qu en - c y [H z] /D am pi ng []) p1 58.1 0.27 90.9 0.2 517 1 34.5 1 257 1 p2 58.1 0.27 90.9 0.2 138 1 22.1 1 69 1 p3 20.2 1 24.5 1 37.1 1 14 1 18.4 1 z1 668 -1 275 -1 187 -e-2 123 8 -1 128 1 1 z2 56.7 1 91.1 1 187 -e-2 14.9 0.93 39 1 z3 26.4 1 28.5 1 40.9 1 14.9 0.93 27.3 1 dc 533 415 5551 -8.642 366.3 Chapter 6: Closed Loop Identification of Feed Drive System 127 6.5.3 Validation on a Virtual 5-Axis Machine Tool The 5-axis spiral contour toolpath shown in Chapter 4 was interpolated with tangen- tial velocity of 75[mm/sec] and acceleration of 21500[mm/sec ], and previously presented contour error models are utilized in computing the errors. Figure 6.13 shows the true and estimated contour errors. As noted, 5-axis contour errors computed from the true and identified model responses are in very good agreement. The discrepancy between the estimated and true contour errors is observed to be less than 10[micron] validating the practicality of the proposed PSO identification method to be used in 5-axis contouring predictions. 6.6 Conclusion This chapter has presented a rapid closed loop identification strategy for machine tool feed drives. The general feed drive model is first introduced that can model feed Figure 6.13: Predicted and actual tracking Performance of the P-PI controlled virtual 5-axis machine tool (Simulation). Chapter 6: Closed Loop Identification of Feed Drive System 128 drives controlled with variety of feedforward and feedback controllers, and a simple linear LS solution is derived to estimate the drive parameters. In case of pole-zero cancellations in the true system, and lack of excitation in the command signal generated by the CNC interpolator causes the ordinary LS solution to produce models with unsta- ble poles. A Particle Swarm Optimization algorithm, which accommodates stability constraints has been developed to tackle the problem, hence, the stability of the identi- fied drive model is guaranteed. The practicality and effectiveness of the proposed PSO identification scheme have been demonstrated in simulation and experimental results, in which the tracking accuracy of virtual and real systems was successfully predicted. Finally, the method has been utilized in estimating the contouring errors that occur during simultaneous 5-axis contour machining applications. Results have proved that the identification strategy can efficiently be used for modeling and predicting the accuracy of 5-axis machine tools. 129 7. Conclusions 7.1 Conclusions Novel smooth optimal trajectory generation algorithms and a precision contour con- trol technique have been developed for 5-axis machine tools. A smooth feedrate planning technique has been formulated for traveling along spline tool-paths in mini- mum time without compromising the tracking and contouring accuracy of the servo system. A novel motion control law has been developed for 5-axis machines, based on analytical modeling of the contour errors that occur during 5-axis motion, which mini- mizes the contouring errors directly instead of axis tracking and capable of achieving high contouring accuracy at high feed speeds. The contributions are summarized as follows: The developed feed planning algorithm has shown that the machining time can be significantly reduced (15-20 [%] shorter cycle times) by optimizing the feed along the tool-path while respecting velocity, acceleration and jerk limits of all five drives. In addition, avoiding the violation of drive saturation limits allows the control law to operate in a linear region, improving the accuracy and performance of the CNC system. Feedrate is optimally planned by overlapping windows with reduced number of variables in real-time on [15]. New analytical contour error estimation models have been developed for 5-axis ma- chine tools in this thesis. Two types of contour errors have been introduced to the literature, the tool tip and orientation contour errors, which are detrimental to part tolerances in 5-axis machining. Proposed contour error models utilize solely the kine- matics of the reference tool-path and actual position of the 5-axis drives, and it is proven Chapter 7: Conclusions 130 in contour following tests that proposed analytical models can estimate contour errors precisely during simultaneous 5-axis machining. Combined with the kinematic model of the 5-axis machine, it has also been shown that all the five axis drives influence tool tip contouring accuracy while the tool axis orientation is controlled only by the two rotary drives. Finally, the models have been utilized to verify that even if tracking accuracies of the all five axis drives are individually greatly improved using high bandwidth and robust sliding mode controllers, the overall contouring error of 5-axis tool paths has still been significant due to the path curvature and the kinematics of the machine [26]. Instead of controlling each drive individually as in the current practice, a multi-input multi-output Sliding Mode Controller that simultaneously controls all the five drives with the objective of achieving high path tracking accuracy has been introduced and experimentally verified in 5-axis contouring tests. The proposed Sliding Mode Control- ler consists of two MIMO SMCs working side by side: The first SMC computes the tool tip contour error using the developed model and generates coupled control action to the Cartesian drives in order to minimize the deviation of the tool tip from the desired trajectory. The second SMC, in parallel, computes the orientation contour error and minimizes it by coupling only the rotary drives. The complete coupled 5-axis sliding mode contour controller has been constructed by proving the stability of the controllers in joint operation. As a result, during high speed 5-axis machining tool tip contour errors have been reduced by 20-30 % using the proposed controller, as compared to the indi- vidual axis based sliding mode controlled case [26],[43]. A closed loop identification method for CNC drive systems has been developed, which can be used in 5-axis contouring simulations when combined with the models Chapter 7: Conclusions 131 developed in the thesis. The technique is based on exciting the feed drive using a series of short motions, and the reference and actual position, velocity and acceleration states of the drive are measured. A general model of the closed loop feed drive system that can represent servo dynamics controlled with various controllers has been fitted to the identification data. In order to avoid identification of an unstable model, the identifica- tion task has been cast as a constrained minimization problem in which the frequency and damping ratio bounds are imposed on the closed loop pole locations, and the nonli- near minimization problem has been solved using a Particle Swarm Optimization algorithm. Simulation and experimental results have shown that drive models con- structed with closed loop identification can be used to predict the contouring accuracy of 5-axis machine tools in a virtual environment with sufficient accuracy. 7.2 Future Research Directions The trajectory generation and control laws developed in this thesis are directly ap- plicable to 5-axis machine tools with different kinematic configurations such as table, spindle tilting or hybrid. Extension of the algorithms to more complex kinematic confi- guration such as parallel kinematic machine tools or robots, is possible, and will allow their benefits to be exploited in a larger variety of applications, such as high speed machining of dies and molds or laser cutting and welding operations. The feedrate scheduling algorithm increases the productivity especially for rough milling or semi-finishing where rapid feed alterations are not vital to part tolerances. In finish milling, on the other hand, the tangential feedrate must be scheduled smoothly where the fluctuations must be kept minimal to avoid leaving tool deflection marks on the part surface. One practical solution would be to consider the tool-path section by Chapter 7: Conclusions 132 section, and the maximum local feed for each section has to be identified. Constant feed sections can be connected to each other with smooth acceleration transients. This will result in a suboptimal solution increasing the cycle time 10~15[%] as compared to the currently developed method, but it will in return reduce the computational cost of planning and allow wider practical application to machine tools. On the other hand, the key to extending the developed time optimal feed profiling to robots and hexapods relies on incorporating drive constraints correctly into the algorithm through proper modeling of their non-Cartesian kino-dynamic chain. The contour control algorithms on the other hand, improve the contouring accuracy by two-fold on a 5-axis machine tool. In order to further increase the stability and precision of the developed algorithms the flexible and dynamic coupling of the drive system has to be considered, at least in the stability analysis and the performance evalua- tion in the control law. A simpler and likely more practical approach may be to use the developed contour error models for a given tool-path and estimate the contour error in a look ahead fashion. In this way, the reference position commands to the drives can then be compensated in order to reduce the contour errors. In parallel, the tool deflection is another issue that contributes to large contour errors in 5-axis machining. Simplified process models are required to be used to estimate the deflection of the tool in real time, and its contribution to contour errors can then be compensated for either in the axis controllers, or as mentioned, in the reference trajectory. The extension of the control algorithms to robots and hexapods is vague through incorporation of their kinematics and dynamics, which will generate further research opportunities to address the new problems. 133 Bibliography [1] Bobrow, J. E., Dubowsky, S., and Gibson, J., 1985, “Time-Optimal Control of Robotic Manipulators Along Specified Paths,” Int. J. Robot. Res., 4, pp. 3–17. [2] Shin, K., McKay, N., 1985, “Minimum-Time Control of Robotic Manipulators with Geometric Path Constraints,” IEEE Trans. Autom. Control, 30, pp. 531–541. [3] Shiller, Z., Lu, H. H., 1990, “Robust Computation of Path Constrained Time Optimal Motions,” IEEE International Conference on Robotics and Automation, Cincinnati, OH, May 13–18, IEEE Comput. Soc. Press, Los Alamitos, CA, Vol. 1, pp. 144–149. [4] Shiller, Z., Lu, H. H., 1992, “Computation of Path Constrained Time Optimal Motions with Dynamic Singularities,” ASME J. Dyn. Syst., Meas., Control, 114, pp. 34–40. [5] Shin, K., McKay, N., 1986, “A Dynamic Programming Approach to Trajectory Planning of Robotic Manipulators,” IEEE Trans. Autom. Control, 31, pp. 491– 500. [6] Gourdeau, R., Schwartz, H. M., 1989, “Optimal Control of a Robot Manipulator Using a Weighted Time-Energy Cost Function,” IEEE Conference on Decision and Control, Tampa, FL, December 13–15, IEEE, Piscataway, NJ, Vol. 2, pp. 1628–1631. [7] Smith, T. S., Farouki, R. T., Timar, S. D., and Boyadjieff, G. L., 2005, “Algo- rithms for Time-Optimal Control of CNC Machines Along Curved Tool-paths,” Rob. Comput.-Integr. Manufact., 21, pp. 37–53. 134 [8] Renton, D., and Elbestawi, M. A., 2000, “High Speed Servo Control of Multi-axis Machine Tools,” Int. J. Mach. Tools Manuf., 40, pp. 539–559. [9] Dong, J., and Stori, J. A., 2006, “A Generalized Time-Optimal Bidirectional Scan Algorithm for Constrained Feedrate Optimization” Journal of Dynamic Systems, Measurement and Control, Transactions of the ASME, Vol. 128, No. 2, pp. 379- 390. [10] Tarkiainen, M., and Shiller, Z., 1993, “Time Optimal Motions of Manipulators with Actuator Dynamics,” IEEE International Conference on Robotics and Auto- mation, Atlanta, GA, May 2–6, IEEE Comput. Soc. Press, Los Alamitos, CA, Vol. 2, pp. 725–730. [11] Piazzi A., and Visioli A., 1998, “Global Minimum-Time Trajectory Planning of Mechanical Manipulators Using Interval Analysis”, Intl. J. of Control, Vol. 71, No. 4, pp. 631-652. [12] Constantinescu, D., and Croft, E. A., 2000, “Smooth and Time Optimal Trajecto- ry Planning for Industrial Manipulators Along Specified Paths,” J. Rob. Syst., Vol. 17, pp. 233–249. [13] Weck, M., Meylahn, A., and Hardebusch, C., 1999, “Innovative Algorithms for Spline-Based CNC Controller”, Prod. Eng. Research and Development in Germa- ny; Annals of the German Academic Society for Production Engineering, Vol. VI/1, pp. 83-86. [14] Altintas, Y., and Erkorkmaz, K., 2003, “Feedrate Optimization for Spline Interpo- lation in High Speed Machine Tools,” CIRP Annals, Vol. 52, pp. 297–302. 135 [15] Sencer, B., Altintas, Y., Croft, E.A., 2007, “Feed Optimization for Five-Axis CNC Machine Tools with Drive Constraints”, Int. Journal of Machine Tool and Manufacture, Vol., 48, Iss., 7-8, pp., 733-745. [16] Koren, Y., 1983, Computer Control of Manufacturing Systems, McGraw-Hill, New York. [17] Koren, Y., Lo, C. C., 1980, “Cross-Coupled Biaxial Computer Control for Manu- facturing Systems”, Journal of Dynamic Systems, Measurement, and Control, Vol. 102, pp. 265-272. [18] Koren, Y., Lo, C.C., 1991, "Variable-Gain Cross-Coupling Controller for Con- touring", Annals of CIRP, Vol. 40/1/1991, pp. 371-374. [19] Erkorkmaz, K., Altintas, Y., 1998, High Speed Contouring Control Algorithm for CNC Machine Tools, Proceedings of the ASME Dynamic Systems and Control Division, ASME International Mechanical Engineering Congress and Exposition, DSC 64, pp. 463–469. [20] Chiu, G. T. C., Tomizuka, M., 2001, “Contouring Control of Machine Tool Feed Drive Systems: A task coordinate frame approach”, IEEE Transactions on Control Systems Technology, Vol. 9, No.1, pp.130-139. [21] Yeh, S.-S., Hsu, P. L., 2002, “Estimation of the Contouring Error Vector for the Cross-coupled Control Design”, IEEE/ASME Transactions on Mechatronics, Ver. 7, No. 1, pp. 44-51. [22] Peng, C. C., Chen, C. L., 2007, “Biaxial Contouring Control with Friction Dy- namics using a Contour Index Approach”, International Journal of Machine Tools and Manufacture, Vol. 47, No.10, pp.1542-1555. 136 [23] Katz, R., Yook J., Koren. Y., 2004, “Control of a non-Orthogonal Reconfigurable Machine Tool”, Journal of Dynamic Systems, Measurement and Control, Trans- actions of the ASME, Vol. 126, No.2, pp.397-405. [24] Lin, R. S., Hu, C. W., 1999, “Modeling of the Dynamic Contour Error of Five- Axis CNC machine tools”, ASME Dynamic Systems Control and Control Divi- sion, Ver. 67, pp. 861-868. [25] Lo, C. C., 2002, “A Tool-path Control Scheme for Five-Axis Machine Tools”, International Journal of Machine Tools and Manufacture, Vol.42, No.1, pp.79-88. [26] Sencer, B., Altintas, Y., Croft, E. A., “Modeling and Control of Contouring Errors for Five-Axis Machine Tools. Part I – Modeling”, Transactions of ASME Journal of Manufacturing Science and Engineering, Vol. 131, Iss. 3, 8 pages. [27] Koren, Y., 1991, “Advanced Controllers for Feed Drives”, CIRP Annals, Vol. 41, No. 2, pp. 689-698. [28] Astrom, K., Wittenmark, B., Computer-Controlled Systems: Theory and Design (3rd ed.), Prentice-Hall, Englewood Cliffs, NJ, 1997. [29] Tomizuka, M., 1987, “Zero Phase Error Tracking Algorithm for Digital Control”, Transactions of ASME Journal of Dynamic Systems, Measurement, and Control, Vol.109, pp.65–68. [30] Tung, E.D., Tomizuka, M., 1993, “Feedforward Tracking Controller Design Based on the Identification of Low Frequency Dynamics”, Transactions of ASME Journal of Dynamic Systems, Measurement and Control, Vol.115, No.3, pp.348- 356. 137 [31] Tsao, T. C., Tomizuka, M., 1987, “Adaptive Zero Phase Error Tracking Algo- rithm for Digital Control”, Transactions of ASME Journal of Dynamic Systems, Measurement and Control, Vol. 109, No. 4, pp. 349-354. [32] Utkin, V. I., 1977, “Variable Structure Systems with Sliding Modes”, IEEE Transactions on Automatic Control, Vol.22, No.2, pp.212-222. [33] Slotine, J. J. E., Li, W., 1988, “Adaptive Manipulator Control: A Cast Study”, IEEE Transactions on Automatic Control, Vol.33, No.11, pp.995-1003. [34] Stepanenko, Y., Cao, Y., Su, C. Y., 1998, “Variable structure control of robotic manipulator with PID sliding surfaces”, International Journal of Robust and Non- linear Control, Vol.8, No.1, pp.79-90. [35] Altintas, Y., Erkorkmaz, K., Zhu, W. H., 2000, “Sliding Mode Controller Design for High Speed Feed Drives”, Annals of CIRP, Vol.49, No.1, pp.265-270. [36] Rodriguez-Angeles, A., Nijmeijer, H., 2005, “Synchronizing Tracking Control for Flexible Joint Robots via Estimated State Feedback”, Journal of Dynamic Sys- tems, Measurement and Control, Transactions of the ASME, Vol.126, No.1, pp.162-172. [37] Su, Y., Sun, D., Ren, L., Mills, J. K., 2006, “Integration of Saturated PI Syn- chronous Control and PD Feedback for Control of Parallel Manipulators”, IEEE Transactions on Robotics, Vol.22, No.1, pp.202-207. [38] Ren, L., Mills, J. K., Sun, D., 2006, “Adaptive Synchronized Control for a Planar Parallel Manipulator: Theory and experiments”, Transactions of ASME Journal of Dynamic Systems, Measurement and Control, Vol.128, No.4, pp.976-979. 138 [39] Erkorkmaz K., Altintas Y., 1998, High Speed Contouring Control Algorithm for CNC Machine Tools, Proceedings of the ASME Dynamic Systems and Control Division, ASME International Mechanical Engineering Congress and Exposition, DSC 64, pp.463–469. [40] Chiu, G.C., Tomizuka, M., 2001, “Contouring Control of Machine Tool Feed Drive Systems: A Task Coordinate Frame Approach”, IEEE Transactions on Con- trol Systems Technology, Vol.9, No.1, pp.130-139. [41] Peng, C. C., Chen, C. L., 2007, “Biaxial Contouring Control with Friction Dy- namics Using a Contour Index Approach”, International Journal of Machine Tools and Manufacture, Vol.47, No.10, pp.1542-1555. [42] Chen, S. L., Liu, H. L., Ting, S. C., 2002, “Contouring Control of Biaxial Systems Based on Polar Coordinates”, IEEE/ASME Transactions on Mechatronics, Vol.7, No.3, pp.329-345. [43] Sencer, B., Altintas, Y., “Modeling and Control of Contouring Errors for Five- Axis Machine Tools. Part II – Precision Contour Controller Design”, Transactions of ASME Journal of Manufacturing Science and Engineering, Vol. 131, Iss. 3, 10 pages. [44] Koren, Y., 1983, Computer Control of Manufacturing Systems, McGraw-Hill, NY. [45] Armstrong, H.B., Dupont, P., Canudas De Wit, C., 1994, “A Survey of Models, Analysis Tools and Compensation Methods for the Control of Machines with Friction,” Automatica, Vol. 30, No. 7, pp. 1083-1138. 139 [46] Lee, H. S., Tomizuka, M., 1996, “Robust Motion Controller Design for High- Accuracy Positioning Systems,” IEEE Transactions on Industrial Electronics, Vol. 43, No. 1, pp. 48-55. [47] Erkorkmaz, K., Altintas, Y., 2001, “High Speed CNC System Design: Part II - Modeling and Identification of Feed Drives,” International Journal of Machine Tools and Manufacture, Vol. 41, No. 10, pp. 1487-1509. [48] Varanasi, K.K., Nayfeh, S.A., 2004, “Dynamics of Lead-Screw Drives: Low- Order Modeling and Experiments,” Transactions of ASME, Journal of Dynamic Systems, Measurement, and Control, Vol. 126, No. 2, pp. 388-396. [49] Erkorkmaz, K., Kamalzadeh, A., 2006, “High Bandwidth Control of Ball Screw Drives,” Annals of CIRP, Vol. 55, No. 1, pp. 393-398. [50] Kao, J. Y., Yeh, Z. M., Tang, Y.S., Lin Y.S., 1996, “A Study of Backlash on the Motion Accuracy of CNC Lathes,” International Journal of Machine Tools and Manufacture, Vol. 36, No. 5, pp. 539–550. [51] Cuttino, J.F., Dow, T.A., Knight, B.F., 1997, “Analytical and Experimental Identification of Nonlinearities in a Single-Nut Preloaded Ball Screw”, Transac- tions of ASME Journal of Mechanical Design, Vol. 119, No. 1, pp. 15-19. [52] Donmez, M. A., Blomquist, D. S., Hocken, R. J., Liu C. R., Barash, M. M., 1986, “A General Methodology for Machine Tool Accuracy Enhancement by Error Compensation,” Precision Engineering, Vol. 8, No. 4, pp. 187-196. [53] Ljung, L., 1998, System Identification: Theory for the User, 2nd Edition, Pren- tice-Hall of Canada Ltd. 140 [54] Wang, L., Cluett, W.R., 1994, “System Identification Based on Closed-loop Step Response Data”, IEE Proc.-Control Theory Appl., Vol. 141, No. 2, pp. 107-110. [55] Silva, R., Barbaro, D., Leon de la Barra, B.A., 2006, “Closed-Loop Process Identification under PI Control: A Time Domain Approach, Ind. Eng. Chem. Res. Vol. 45, pp. 4671-4678. [56] Park, H. I., Sung, S. W., Lee, I.B., 1997, “On-line Process Identification Using the Laguerre Series for Automation Tuning of the Proportional-Integral- Derivative Controller”, Ind. Eng. Chem. Res. Vol. 36, pp. 101-111. [57] Fliess, M., Sira-Ramirez, H., 2003, “An Algebraic Framework for Linear Identifi- cation”, Control, Optimization and Calculus of Variations, Vol. 9, pp. 151-168. [58] Sira-Ramirez, H., Barrios-Cruz, E., Marquez, R. J., 2007, “Fast Adaptive Trajec- tory Tracking Control of a Completely Uncertain DC Motor via Output Feedback”, 46th IEEE Conference on Decision and Control New Orleans, LA, USA, pp. 4197-4202. [59] Erkorkmaz, K., Wong, W., 2006, “Rapid Identification Technique for Virtual CNC Drives,” Proc. CIRP 2nd Intl. Conf. on High Performance Cutting (HPC'06), Vancouver, BC, June 12-13. [60] Erkorkmaz, K., Wong, W., 2007, “Rapid Identification Technique for Virtual CNC Drives”, International Journal of Machine Tools and Manufacture, Vol. 47, pp. 1381-1392. [61] Wong, W., Erkorkmaz, K., 2007, “Constrained Identification of Virtual Drive Models Using a Genetic Algorithm”, Proc. CIRP Sponsored International Confe- 141 rence on Smart Machining Systems, National Institute of Standards and Technol- ogy (NIST), USA, March 13-15. [62] Eberhart, R. C., Kennedy, J., “A New Optimizer Using Particle Swarm Theory,” in Proceedings of 6th Int. Symp. Micromachine Human Sci., vol. 1, 1995, pp. 39- 43. [63] Devanit, R.S., Hartenberg, J., 1955, “A Kinematic Notation for Lower Pair Mechanisms Based on Matrices”, ASME Journal of Applied Mechanics 22(2) 215-221. [64] Bohez, E. L. J., 2002, “Five-axis Milling Machine Tool Kinematic Chain Design and Analysis”, International Journal of Machine Tools and Manufacture, Vol.42, No.4, pp.505-520. [65] Jung, Y. H., Lee, D. W., Kim J. S., Mok, H. S., 2002, “NC Post-processor for 5- Axis Milling Machine of Table-rotating/tilting Type”, Journal of Materials Processing Technology, Vol. 130-131, pp. 641-646. [66] Tutunea-Fatan, O. R., Feng, H. Y., 2004, “Configuration analysis of Five-Axis Machine Tools Using a Generic Kinematic Model”, International Journal of Ma- chine Tools and Manufacture, Vol.44, No.11, pp.1235-1243. [67] Erkorkmaz, K., Altintas, Y., 2001, “High Speed CNC System Design: Part I – Jerk Limited Trajectory Generation and Quintic Spline Interpolation,” Interna- tional Journal of Machine Tools and Manufacture, Vol. 41, No. 9, pp. 1323-1345. [68] Fleisig, R. V., Spence, A. D., 2001, “Constant Feed and Reduced Angular Accele- ration Interpolation Algorithm for Multi-Axis Machining”, CAD Comp. Aided Design, Vol. 33, No. 1, pp. 1-15. 142 [69] Erkorkmaz, K., Altintas, Y., 2003, “Quintic Spline Interpolation with Minimal Feed Fluctuation,” Proceedings of the ASME Manufacturing Engineering Divi- sion, 2003 ASME International Mechanical Engineering Congress and Exposition, MED-IMECE2003-42428. [70] Do Carmo, M. P., “Differential Geometry of Curves and Surfaces”, Englewood Cliffs, NJ. Prentice-Hall, 1976. [71] Coleman, T., Branch, M. A., and Grace, A., 1999, Optimization Toolbox for Use with MATLAB®, User’s Guide Version 2, The MathWorks Inc. [72] Altintas, Y., Erol, N. A., 1998, “Open Architecture Modular Tool Kit for Motion and Machining Process Control”, Annals of CIRP, Vol. 47, No. 1, pp. 295-300. [73] Slotine, J. J., Li, W., Applied Nonlinear Control, Prentice Hall, NJ, 1991. [74] Nayan, R. S., Amit, K., Swagatam, D., Ajith, A., 2007, “A Closed Loop Stability Analysis and Parameter Selection of the Particle Swarm Optimization Dynamics for Faster Convergence”, 2007 IEEE Congress on Evolutionary Computation , pp.1769-1776. [75] Zhang, Y., Li, S., 2007, “Parameter Identification for Wiener Model Using Particle Swarm Optimization with a Case Study”, Proceedings of the IEEE Inter- national Conference on Automation and Logistics, pp. 1725-1730. [76] Visakan, K., Kirusnapillai, S., Peter, J. F., 2006, “Stability Analysis of the Particle Dynamics in Particle Swarm Optimizer”, IEEE Transactions on Evolutionary Computation, Vol. 10, No. 3, pp. 245-255.