UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Five-axis trajectory generation methods Sencer, Burak 2005

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
831-ubc_2005-0631.pdf [ 7.17MB ]
Metadata
JSON: 831-1.0080755.json
JSON-LD: 831-1.0080755-ld.json
RDF/XML (Pretty): 831-1.0080755-rdf.xml
RDF/JSON: 831-1.0080755-rdf.json
Turtle: 831-1.0080755-turtle.txt
N-Triples: 831-1.0080755-rdf-ntriples.txt
Original Record: 831-1.0080755-source.json
Full Text
831-1.0080755-fulltext.txt
Citation
831-1.0080755.ris

Full Text

FIVE-AXIS TRAJECTORY GENERATION METHODS By Burak Sencer B. Sc., Istanbul Technical University, 2003  A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF APPLIED SCIENCE  in  THE FACULTY OF GRADUATE STUDIES (MECHANICAL ENGINEERING)  THE UNIVERSITY OF BRITISH COLUMBIA September 2005  © Burak Sencer, 2005  Abstract Five-axis machining has been widely adopted manufacturing parts with sculptured surfaces in the aerospace, die and mold and other industries. Five-axis machine tools have a kinematic structure compromising three Cartesian axes couples with two rotary joints enabling the tool to reach any given position and the orientation on the workpiece surface. Five-axis motion allows the orientation of the tool along the toolpaths with varying curvatures and sculptured surfaces. However, optimal trajectory generation methods are required in achieving a smooth coordinated motion of the tool tip and tool axis orientation, which is essential in reducing the fluctuations in the feed, accelerate and jerk of the drives. In this thesis, a toolpath generation algorithm has been developed for smooth trajectories based on compressing the C° continuous Numerical Control (NC) motion commands into continuous b-spline format. With the aid of the b-spline toolpath representation, the discontinuities and unwanted acceleration harmonic in the multi-axis motion are reduces. The tool position and orientation vectors are interpolated continuously enabling a smooth and synchronized motion of the drives. In order to achieve higher smoothness in the change of feed direction, the spline toolpath is optimized by solving the linear quadratic minimization problem during the least squares fit. The resultant angular motion of the tool is further smoothed by minimizing the integral square of third derivative of the orientation spline using a non-linear optimization technique. The Feed Correction Polynomial algorithm is implemented in the real time interpolation of the generated b-spline toolpath yielding a consistent feedrate profile during contouring operations. It is shown that feed correction method exhibits reduced feedrate fluctuations when compared to the widely used Taylor Series approximations. A jerk continuous feedrate profile in terms of C quintic spline is adapted for five-axis 3  trajectory generation. The feedrate profile is generated utilizing the minimum jerk criteria for the smoothness of the motion and tracking accuracy. The feed profile is then optimized to achieve minimum cycle time while adhering to the machine tool's velocity, acceleration and jerk limitations. The developed trajectory generation algorithms have been simulated using the kinematics of a widely used five-axis machine tool promising potential applications in rough machining.  ii  Table of Contents  Abstract  ii  Table of Contents  iii  List of Tables  vi  List of Figures  vii  Nomenclature  ix  Acknowledgment  xiv  Chapter 1 Introduction  1  Chapter 2 Literature Review  3  2.1. Overview  3  2.2. Spline Properties  5  2.3. Spline Fitting for Sculptured Machining  8  2.4. Trajectory Generation for 3 Axis Motions  12  2.5. Critical Issues on 5-Axis Spline Trajectory Generation  13  Chapter 3 Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  18  3.1. Introduction  18  3.2. B-Spline Background  20  3.3. Double B-Spline Polynomial Representation  24  3.4. Generation of Lower B-Spline  27  3.4.1. Tolerance Definition and Parameter Correction 3.5. Generation of the Upper Spline  33 39  iii  3.6. Smoothing the Double Spline Toolpath  42  3.6.1. Introduction  42  3.6.2. Smoothing the Lower B-Spline  43  3.6.3. Angular Jerk Minimization  47  3.7. Simulation Results  50  3.8. Summary  63  Chapter 4 Interpolation and Feedrate Profiling for 5-Axis Machining  64  4.1. Introduction  64  4.2. Interpolation with Minimal Feed Fluctuation  64  4.2.1. Taylor Series Approximation  66  4.2.2. Feed Correction Polynomial Method  70  4.2.3. Iterative Solution of the Spline Parameter  77  4.3. Feedrate Profiling for 5-Axis Machine Tools  79  4.3.1. Minimum Jerk Feed Profile  81  4.3.2. Velocity, Acceleration and Jerk Constraints  84  4.3.3. Inverse Kinematics of a Rotary Tilting Table 5-Axis Machine Tool  84  4.3.3.1 Rotational Motions  88  4.3.3.2 Cartesian Motions  91  4.3.4. Formulation of the Constraints  95  4.3.5. Solution of the Optimization Problem  97  4.4. Simulation Results  99  4.4.1. Comparison of the Interpolation Methods  99  4.4.2. Comparison of Constant and Optimized Feedrate Profiles  102  4.4.3. Summary  105  Bibliography  111  Appendix A Formulation of the Minimum Jerk Feed Profile  120  Appendix B Formulation of Gradients of Feedrate Optimization  126  B.l. Gradient of the Objective  125 iv  B.2. Gradient of the Feed Profile Coefficients B.3. Gradient of the Kinematic Profiles B.4. Gradients of the Constraints  List of Tables Table 3.1: Comparison of Curvatures  56  Table 4.1: Rotary Axis Angles  91  Table 4.2: Summary of Interpolation Results  100  Table 4.3: Kinematic Limits of 5-Axis Machine Tool  103  vi  List of Figures Figure 2.1: Reference Trajectory Generation  3  Figure 2.2: Quintic Polynomial Spline  6  Figure 2.3: Bezier Curve  7  Figure 2.4: B-spline curve  8  Figure 2.5: Tool Tip and Tool Orientation Vectors  14  Figure 2.6: a) Position Spline — b) Orientation Spline  15  Figure 3.1: Flow Chart of Optimum B-spline Toolpath Generation  19  Figure 3.2: B-spline Toolpath  21  Figure 3.3: B-spline Basis Functions  23  Figure 3.4: Double B-Spline Toolpath Representation  25  Figure 3.5: Upper and Lower B-splines  26  Figure 3.6: Initial Spline Parameter Calculation  28  Figure 3.7: Tolerance Definition for Lower B-spline Generation  34  Figure 3.8: Spline Parameter Correction  35  Figure 3.9: Zero Tolerance Definitions a) Zero Cosine — b) Point Coincidence  36  Figure 3.10: Lower B-spline Fitting Algorithm  38  Figure 3.12: Double B-spline Generation Algorithm  41  Figure 3.13: Tool Axis Vector  47  Figure 3.14: Test Double B-spline Toolpaths Generated from Discrete Points  52  Figure 3.15: Fitting Errors of the B-Spline Toolpaths  53  Figure 3.16: Multi Objective Least Squares Smoothening of a B-spline - a) Regular Least Squares Fit, b-d) Optimized Least Squares Fit Comparisons  55  Table 3.1: Comparison of Curvatures  56  Figure 3.17: 2 and 3 Derivative Comparison of Smoothened B-splines  56  Figure 3.18: Comparison of Optimized Lower B-Spline  58  Figure 3.19: Angular Kinematic Profiles  60  Figure 3.20: Rotary Axis Kinematic Profiles  61  Figure 3.21: Linear Axis Kinematic Profiles  62  Figure 4.1: Real-time Spline Interpolation  65  nd  rd  vii  Figure 4.2: Basis Function Computation a) Recursive Calculation of Basis Functions — b) Recursive Calculation of Basis Function Derivatives  69  Figure 4.3: B-spline Knot Spans and Segmentation  71  Figure 4.4: Feed Correction Polynomial for B-spline Toolpath  73  Figure 4.5: Real Time Interpolation Methods  77  Figure 4.6: Overview of Feedrate Profiling  80  Figure 4.7: Kinematic Profiles on the B-spline Toolpath  82  Figure 4.8: Kinematic Chain of the Rotary Tilting Table 5-Axis Machine Tool  85  Figure 4.9: a) Translation Motion — b) Rotational Motion  86  Figure 4.10: a) Table Orientation Before Cutting — b) Table Orientation During Cutting 89 Table 4.1 Rotary Axis Joint Angles  91  Figure 4.11: Kinematic Transformation  94  Figure 4.12: Constraint Evaluation  98  Figure 4.13: B-spline Toolpath  101  Figure 4.14: Feedrate Fluctuations  101  Figure 4.15: B-spline Toolpath  103  Figure 4.16: Optimized and Unoptimized Feedrate Profiles  104  Figure 4.17: Optimized and Unoptimized Cartesian Axis Profiles  104  Figure 4.18: Optimized and Unoptimized Rotary Axis Profiles  105  viii  Nomenclature  A = A{u)  = tool axis spline  A(,B[,...H[  = feed correction polynomial coeff.  4M--M  = normalized feed correction polynomial coeff.  A>B  = quintic spline displacement profile coeff.  c  k  k  a ,b ,...,c  = quintic spline displacement profile normalized coeff.  A ,A ,A  = x,y,z axis components of tool axis spline  k  k  x  y  k  2  X max ' Ymax> Z a  a  ermax max •e>max> "6. max ' "6 max > >y>-Z, a  a  =  C  x  K C axis acceleration limits  C  = vector of inequality constraints  D  = vector between estimated and actual data points  e  = error, position and parameter  Cmax  = maximum tool tip positioning error  ee  = orientation error  &6max  = maximum tool orientation error  fmax  = maximum tangential feedrate  H  = relative distance between lower and upper B-splines  TjRef 11  QT  reference frames  = homogeneous transformation matrices between tool and  H  Ref  WP  11  and reference frames  =homogeneous transformation matrices between workpiece  JK JL JL  Jy K-~x,y,z  = objective function for smoothing = objective function of least squares fit = augmented objective function for least squares = objective function for minimum jerk = axis components of smoothing objective function = axis components of the least squares objective function  J  = objective function for minimum angular jerk  0  J  = objective function for minimum jerk disp. profile  J  — augmented objective function for minimum jerk disp profile  J  = segment objective function for minimum Jerk disp profile  v  v  v  >Y> > A, C aXeS jerk HmitS Z  JXmax' JYmax' JZ max' J6A max • J8C max ~ > X  K  = weight matrix  L  = constraint matrix  M  = number of fitting data points - 1  m  = number of knots - 1  M(t), M(t), M(t)  = x, ,y, z, A, C axes velocity, acceleration and jerk profiles  M M, M  = geometric derivatives of the toolpath  n  = number of control points -1  s  ss  sss  N ,N ,N ,...,N  = basis functions  Ns  = number of segments (internal knot spans)  j0  n  i2  ip  N ,N ,N  = 1 ,2  O - 0(u)  = Orientation spline  O  = reference orientation  u  m  st  mu  r  nd  ,3  rd  derivatives of the basis functions  O ,O , 0[  = x, y, z axis components of reference orientation  O' '  = fixed tool orientation  r  r  x  y  00  OT OT Ot  = tool frame axial offset components  O, O , O  = 1 , 2 and 3 derivative of the orientation Spline  O , O ,O  = x,y,z axis components of orientation spline  Pj  = control points of lower(tool-tip) B-spline  p  = degree of the B-spline  p{u)  = Tool tip (Lower) B-spline  x  u  x  y  z  uu  st  y  uuu  z  u >uu'uuu  p  p  P ,P ,P  p  x  Z  nd  =  1  rd  2  nd a n d 3 f d  derivative of Lower B-spline  = control points of lower B-spline x,y,z axis components  Qi  = control points of upper B-spline  Q(u)  = upper B-spline  q'o, q[ , • • •, q[, q'  = fitting data points  M  q'  = x axis component of estimated data points  q'x > q'y„ > q'z k  Q '  = k  Quu • Quuu  u  =  =  Qx'QyQz  r  fi > rst  X  >Y'  Z a x  i components of actual Data Points  s e c o n a l  s  a n  d third derivative of Upper B-spline  control points of upper B-spline x,y ,z axis components  = total number of constraints  Re f ,Ref,  Re f  x  = reference frame axial offset components  z  R ,R ,R  = rotary ( C ) frame axial offset components  5 = s(t)  = displacement profile  s (v)  = velocity profile  s (v)  = tangential acceleration profile  s' (v)  = jerk profile  x  y  z  k  k  k  s  fmai'  $ final > 'fmai s  =  f " l velocity, acceleration and jerk conditions ma  s , s , s]  = initial velocity, acceleration and jerk conditions  S  = segment arc travel lengths  injt  ln!t  nil  k  Sp Sp Sp  = spindle frame axial offsets  S  = total arc length spline  x  y  z  lolal  S S S x.  y.  = table frame axial offsets  2  T ,..., T ,..., T _, 0  k  N  = segment time durations  T  = maximum expected segment time duration  Ts  = servo control loop sampling interval  T  - total travel time  max  a  T >T ,T  = tilting (A) frame axial offsets components  U  = knot vector of the B-spline  x  y  z  xi  u  = spline parameter  u  = parameter values corresponding to data points  u  = estimated B-spline parameter  u  = previous spline parameter  u  - tested spline parameter  V = v{f)  = tangential feedrate  prev  lesl  x ax' Y x> z x' e max> e ax=  v  v  v  v  m  ma  ma  v  A  cm  a x i s  velocity limits  Wp Wp Wp  = workpiece frame axial offsets  X,Y,Z  = x,y,z axis position commands or indices  x  y  z  x ,y ,z  = x,y,z axis components of Upper B-spline  u  = previous position commands  x', y', z'  = x,y,z axis components of Lower(tool-tip) B-spline  h  h  h  prev  x ,y , z , tesl  lest  les  = tested position commands  xii  = regressor matrix  0  = constraint vector for quadratic programming A  = vector of Lagrange multipliers  £  = relaxation tolerance for nonlinear optimization = measure of the Euclidean distance = zero cosine value  s  2  o  = vector of upper B-spline optimization variables  a  = fitting effect coefficient  r  = smoothing effect coefficient  £  = relaxation tolerance  a  = normalized arc length travel of segment  h  = vector of Feed Correction Poly, optimization variables  K  = normalized time duration  Of  = vector of feed correction polynomial coefficients = minimum jerk disp. profile optimization variables  o  = minimum jerk orientaiton profile optimization variables  X  = spline displacement profile time parameter  V  = normalized time parameter  h  o  A  = tilting axis (A) rotation angle = rotary axis (C) rotation angle  As  = arc increment  Au  = parameter increment  M  = knot  A  = inserted knot  xiii  Acknowledgment I would like to express my sincere appreciation to my research supervisor Dr. Yusuf Altintas for his valuable instruction, guidance, and support which he has provided throughout my graduate studies at University of British Columbia. I would like to thank Nesrin Altintas for her hospitality. I wish to thank my colleagues in the Manufacturing Automation Laboratory for sharing their knowledge and experience with me. I also feel very luckily to gain such an exceptional friend, Doruk. He was with me especially during the tough times and made my life colorful and enjoyable. Finally, I am deeply grateful to my beloved family, my mother Serap and my father Necip for their love, encouragement unweaving support and patience they have provided me throughout my entire life and abroad studies. This thesis and all of my previous success are dedicated to them.  xiv  1  Chapter 1 Introduction The Computer Aided Design (CAD) technology has led to smooth spline representation of space curves and complex free-form surfaces with steep curvatures, commonly used in designing dies, mold and aerospace parts with sculptures surfaces. Bezier splines (b-splines) and non-uniform rational b-splines (NURBS) are extensively used because of their dexterity in modifying and representing curved geometries with high smoothness. Sculptured machining of free-form surfaces by five-axis computer-numerically-controlled (CNC) machine tools has received significant attention in recent years since it provides optimal tool posture and orientation, not possible with three-axis Cartesian machines. In particular, surfaces with varying curvatures are machined on five-axis machining centers where the tool orientation is adjusted with varying curvatures are machined on five-axis machining centers where the tool orientation is adjusted continuously resulting in a smoother surface finish. Furthermore, on such machines, machining of deep cavities and inclined drilling are performed using shorter cutters, leading to less vibrations and higher material removal rates. However, current practice in five-Axis machining approximates the contours with straight line segments yielding large C° continuous NC toolpaths, composed of tool position and orientation vectors. The linear segmentation technique leads to a rough surface finish and forces the machine tool to decelerate, accelerate or stop at each segment junction. In particular, for machines with rotary axes located on the table, the axes have to move and support significant masses during machining. As a result, the acceleration limits of such fiveaxis machines are reduces, increasing the machining cycle times. Tool motion with abrupt deceleration and accelerations cause contouring errors found at locations where the transients such as tangency or curvature discontinuities occur in the reference trajectory. Such discontinuous feed motions may excite the natural modes of the mechanical structure of the mechanical structure and servo control system, characterized by vibration marks on the part surface and degradation of the axis tracking performance. The discontinuous toolpath must, therefore, be replaced with spline representations for smooth trajectory generation spline toolpaths facilitate the delivery of the desired, smooth, high  Chapter 1. Introduction  2  speed feed motion by enforcing acceleration and jerk continuity in the trajectory. In this thesis, the smooth and continuous toolpath requirements are implemented using b-splines that are capable of interpolating tool position and orientation vectors, and generated by compressing the original NC codes. Unfortunately, spline trajectories exhibit undesirable feedrate fluctuations during real time interpolation. When the curves are interpolated at constant parameter increment, fluctuations in the tangential feedrate occur, which reduce the smoothness of the coordinated multi-axis motion. The appropriate spline parameter must be calculated in order to realize the desired feed profile during machining operation, b-spline toolpaths are used with the Feed Correction Polynomial method in order to smoothly parameterize the curved paths for fast and accurate real time interpolation. The feedrate profile of the trajectory must also maximize smoothness for tracking performance and minimize cycle time for productivity. The optimized smoothness can be achieved by minimizing the jerk discontinuity of the overall trajectory. However, due to the kinematic limits of the drives and tracking performance of the servo system, the tangential feedrate profile is generated conservatively, compromising the minimum cycle time objective. The inverse kinematics transformation of the five-axis machine tool structure further diminishes the linear relationship between the tool motion in workpiece coordinates and axis kinematic profiles. As a result, the optimum scheduling of the feedrate can only be achieved by considering the curvature of the toolpath, the change in the orientation and the kinematic transformation, while adhering the machine tool's velocity, acceleration and jerk limitations. By adapting the feedrate profiling to accommodate these non-linear factors, the high speed tracking and contouring accuracy during the execution of the optimal feed motion can be preserved. The thesis is organized as follows: The review of the related literature is presented in Chapter 2, followed by smooth spline toolpath generation method in Chapter 3. The real time interpolation of the b-spline toolpath and the feedrate scheduling algorithm for five-axis Machining is given in Chapter 4. The conclusions and possible future research directions are discussed in Chapter 5.  3  Chapter 2 Literature Review  2.1. Overview  Reference trajectory generation plays a significant role in five-axis computer numerical controlled (CNC) machine tools to acquire the desired part quality. Generation of the reference motion commands to the machine tool's axis drives consist of three main stages: Planning the machining toolpath on the part, imposing the tangential feedrate profile and real time interpolation of the position commands within the servo loop control frequency.  TOOLPATH G E O M E T R Y Toolpath  o CD Q.  Position - P(u)=x(u)i+y(u)j+z(u)k . . Oriehtatioh-0(u)=Ox(u;/H-Oy(ujy+O (t/JAc| -linear -circular -spline(Quintic,Bezier B-spline...) ;  z  O .O  NC Code  FEED MOTION  c o TD  Planning  INTERPOLATION Inverse Kinematics& Real Time Interpolation]  O W  X = x(t),Y = y(t),Z = z(t) A = a(f),C = c(f)  CD CD CD  Q  vy vy  MACHINE TOOL  Figure 2.1: Reference Trajectory Generation  Chapter 2. Literature Review  4  The attempt in generating a smooth trajectory starts at the toolpath planning step. Small discontinuities in reference toolpaths can result in undesirable high frequency harmonics in reference trajectory that excite the natural modes of the mechanical structure or servo control system. The resultant high frequency components can cause actuator saturations and axis tracking errors causing violations of the machining tolerances on the manufactured parts. Furthermore, employing only linear and circular interpolation techniques to machine complex shapes such as dies, molds, turbine blades, and aerospace parts has serious limitations in productivity since the machine tool must decelerate/accelerate or stop between consecutive G codes [42]. Some amount of work has been done to overcome these problems by developing spline tool path definitions for 3 [78, 79, 80] and five-axis machine tools [1, 29, 55]. Spline toolpaths are parametric curves where the spline parameter does not have an explicit relationship to the total travel arc length. Therefore, ensuring the desired feedrate profile on a spline tool path requires precise calculation of the spline parameter increment during real time interpolation. Feed fluctuations posed by improper calculation of the curve parameter increase the contouring errors. A significant amount of research has been done in real time spline interpolation techniques [20, 42, 50]. On the other hand, the kinematic profiles of individual axes are defined by the feedrate profile of the trajectory. Thus, the feedrate profile must avoid any discontinuity in velocity, acceleration and jerk. To enhance the tracking performance of the servo systems, overall jerk content of the trajectory needs to be minimized and to increase the productivity the feed profile must be scheduled for minimum cycle time [3] while respecting the machine tool's kinematic limits. This is a challenge particularly required to be tacked in five-axis high speed machining centers where the kinematic transformations are significant. The literature reviewed about the trajectory generation techniques for multi axis machining is presented in the following sections. This chapter organized as follows: the spline properties used in tool path definition are introduced in Section 2.2. This is followed by a literature survey of the spline fitting techniques for sculptured surface machining. The research on trajectory generation for 3Axis machining is briefly presented in Section 2.3. A literature review on 5-Axis spline interpolation and motion generation is covered in Section 2.4.  Chapter 2. Literature Review  5  2.2. Spline Properties  Splines are powerful tools that have intensively being used to design smooth and complex geometries such as dies, molds, gas tribune blades and aeronautical parts in the computer aided design (CAD) stage. The most conventional method of presenting spline curves is using implicit equations [27] and parametric functions [58]. However, the majority of CAD systems use the parametric forms of curves instead of implicit since it is practical for part design and control of multi axis machine tools where each axis can be driven individually. Independent from their parametric representation, splines in computer aided design and manufacturing (CAM) systems can be grouped in three main categories:  •  Polynomial splines,  •  Bezier Curves and  •  B-splines.  Cubic and Quintic splines are extensively used curve formats included in the polynomial spline category. Their definition bases on exactly interpolating certain design points (control knots) with piecewise polynomials. The polynomial curve segments are constructed in a way that the position and derivative continuity throughout the composite curve is preserved. A part from the cubic polynomial splines, quintic spline has taken emphasis in the area of CNC tool path generation from the researchers because of their advantageous parameterization [78] for real time interpolation and third order (C ) continuity allowing implementation of 3  jerk continuous trajectories [80]. The quintic spline represented in Fig. 2.2 is composed of N - 1 piecewise curve segments S ...,S interpolating the Preference p o i n t s P P _ \ . It is l  N  l  N  well defined in the polynomial form of  P{u) = Ap + B.U 5  + C,u  3  + Diu + EJU + F, 1  where for two axis cartesian case, S = [S t  xi  Sf  and A = [A  yi  t  The unknown coefficients of quintic splines (A ,B ,...,F ) j  I  J  xi  A ^,...,F y  T  = [F  XI  F J. yi  are found by imposing position  (C°), and higher order continuity boundary conditions between consecutive piecewise curve  Chapter 2. Literature Review  6  segments. However, the affect of coefficients on the spline shape is totally unpredictable, which makes the polynomial splines unsuitable for geometric design.  Figure 2.2: Quintic Polynomial Spline  Bezier curves are being extensively exploited in C A D because it merges the comprehensiveness of the parametrical format with the design functionality. Their definition is clearly presented by Piegl [26] in terms of basis (blending) functions B (u) and control jp  points />=[/>„  P fas, y/  ;=o  where p is the degree of the curve, which also specifies the number of total control points included in the curve definition. The blending functions are given in Bernstein form [26] by:  5,,,W =i\(p-if. - 7 ^ : « ' ( l - ^ r ue[0 1]  Chapter 2. Literature Review  7  As shown in Fig. 2.3, the Bezier curve interpolates only the first and last control points. Its complete shape is determined by the positions of the remaining points. Therefore, the shape can easily be adjusted by modification of any control point making it practical for geometric design.  Figure 2.3: Bezier Curve On the other hand, large number of control points is required to generate complex surfaces resulting higher order Bezier curves. This drawback is defeated with b-splines introduced by de Boor [58]. Similar to Bezier curves, b-splines (See Fig. 2.4) are also defined by control points and basis functions. However, the number of control points and the degree can be selected independently allowing designers to work with lower degree curves. Furthermore, the complete spline is segmented into parameter intervals by a knot vector (U), which enables local modification of the curve by using the nearest control points (local modification property). A p order b-spline curve is defined by, th  where N (u) are the new basis functions evaluated recursively using the de Boor algorithm in  [58].  8  Chapter 2. Literature Review  Figure 2.4: B-spline curve  2.3. Spline Fitting for Sculptured Machining  It was a milestone for sculptured machining to represent complex parts with smooth curves, namely quintic, Bezier and the b-splines. But sending this geometric information with minimal lost to the machine tool in smooth recognizable formats such as splines is vital for high speed machining, since higher and smoother feedrate profiles as well as better tracking performance can be obtained. However, the common practice in generating toolpaths for machining sculptured surfaces on 3 and five-axis machining centers bases on approximating the original surface with lines and sending this cutter location data stacked in a G code to the CNC Executive [62, 63]. Previous researchers Koren and Lo [42], and Yang and Kong [86] have demonstrated the following limitations of this linear interpolation method: Velocity discontinuities at the linear segment junctions lead to higher accelerations, poor surface finish, lower part accuracy and longer machining times. Jouaneh et. al. [35] suggested using arc segments to connect the linear segments and to smooth the corner of the tool path in order to obtain the velocity continuity (C ) at the 1  9  Chapter 2. Literature Review  junctions making it easier to track. Similarly, Hua et. al. [61] used only arcs to discretize the curved contours for 3-Axis high speed machining. Yeung [91] pointed out that quintic splines can be used to connect linear segments, which allows acceleration continuity at the transitions and showed that the tracking performance of the servo system is improved tremendously. Instead of linear and arc interpolation techniques, spline refitting and direct spline trajectory generation methods have been widely adapted in machining sculptured surfaces for HSM. Especially, cubic and quintic splines trajectories are utilized for real time toolpath interpolation because of their advantages in real time interpolation. Wang and Yang [78] introduced C cubic and quintic splines enabling acceleration continuous feed profiles. They 2  pointed out that when splines are interpolated with constant chord length increment, feedrate fluctuations occur due to the discrepancy between the chord length parameter and the traveled arc length increment. They brought up the idea of fitting quintic splines with NearArc Length parameterization resulting in smaller feedrate fluctuations and smoother geometric derivatives. Later, Wang et. al. [79, 80] improved this methodology by preserving the continuity of third derivatives at the consecutive spline segment junctions. This improved technique is called Approximately Arc-Length Parameterization with C and reduced the 3  feedrate fluctuations yielding smoother geometric derivatives. Erkorkmaz and Altintas [23] proposed Optimally Arc Length Parameterized quintic splines. This approach is based on optimizing the spline parameters in such a way that the discrepancy between the spline chord length parameter and traveled arc length increment is minimized, resulting less feedrate fluctuations compared to Near Arc Length parameterization. Furthermore, research has been done in interpolating the quintic spline toolpath in real time when the parameterization is arbitrary. Shipitalni et. al. [67] proposed an iterative interpolation technique using the first order Taylor series approximation in order to determine the spline parameter at constant feedrate. However, this approach was not capable of eliminating the feed fluctuations in case of sharp curvature changes. Later, Lin [50] used second order approximation reporting two order of magnitude reduction in the feed fluctuations compared to the first order. Altintas and Erkorkmaz [21] developed an iterative method requiring the numerical solution of the root of a high order polynomial relating the arc increment to the parameter increment. They also used an additional polynomial to express  Chapter 2. Literature Review  10  the nonlinearity between the spline parameter and the arc length, which is called Feed Correction Polynomial [19] . By the aid of this polynomial, the spline parameter can be computed efficiently in real time acquiring minimum feed fluctuation. Besides the quintic polynomial splines, b-splines are widely used in generating parametric surfaces from discrete measurement data in reverse engineering [43, 65, 66], and b-spline trajectories for machine tools [5, 33, 47, 48] and robots [24, 92]. Since b-splines are capable of representing geometries with higher smoothness and continuity [73] compared to quintic polynomials, they are more suitable for high speed machining purposes. Lai et. al. [46] utilized b-splines to generate sculpture surfaces from laser-scan measurement data and extracted toolpaths for 3-Axis machining. He used a non-linear least squares technique to fit the axis splines within acceptable fitting tolerances. Later, Elkeran and El-Baz [17] used the idea to directly generate non uniform rational b-splines (NURBS) linear G code blocks using least squares. They fed the NC unit with NURBS and showed that consistent feed rates can be achieved as compared to the linear interpolation. Yau and Kuo [89] used the same method to compress G codes into b-spline formats, which can be accepted by different machine tool NC units such as Fanuc and Siemens. Even those sophisticated CNC executives are equipped nowadays with build - in functions that compress numerous blocks of G codes into smooth axis b-splines for HSM [70]. However, the quality of the curves as well as their real time interpolation technique is very important for the tracking performance and surface quality. These problems are addressed from two directions: generation of smooth axis splines and accurate interpolation of the b-splines in real time. Generally speaking, smooth curves with less change in the curvature and jerk content are desired since they are easier to track by the servo system especially when traveled at high feedrates. The problem of curve quality is first addressed by Horn [32]. He showed that maximum smoothness of the curve is achieved when the square integral of the curvature is minimized. Hohenberger and Reudig [31] used the curvature cost function and proposed an unconstrained iterative optimization scheme to approximate measurement data with smooth b-spline surfaces. They pointed out that smoother curves can be obtained while sacrificing from the fitting accuracy. Bercover [2, 30] used various energy functionals for curve smoothness based on optimizing the b-spline control points as well as the knot vectors. On the other hand, Vassilev [75] and others [26, 56] studied the b-spline  Chapter 2. Literature Review  11  fitting while pursuing less computationally expensive methods for computer graphics. They investigated the affect of cost functions formulated by square sum of derivatives with respect to the curve parameter for smoothness. The oscillations on the b-spline are eliminated by adding energy cost functions including first and second derivatives. In these methods, the curve smoothness and fitting accuracy is controlled using several parameters, which are tuned by trial and error during approximation making it unfavorable to implement. The real time interpolation of the b-splines has become an important research area in relation to curve fitting. Greenway and Zhang [92] first put in practice a b-spline curve motion interpolator for a six-axis robot based on Taylor's first-order expansion and compared the feed consistency to linear interpolators. Yang and Park [87] implemented a b-spline interpolator for wire electrical discharge machining (WEDM) also based on Taylor's first order expansion on a open architecture control system. Zhiming et. al. [93] used the first and second order Taylor's expansion with variable feedrate for real time b-spline interpolation. They adjusted the feedrate with respect to the curvature of the toolpath in order to consider the machine tool's dynamics. Tsai et. al. [74] compared several real time interpolation methods of b-splines such as the first and second order Taylor expansion, fourth order Runge - Kutta and an iterative parameter correction scheme. They indicated that second order • Taylor's expansion gives less feedrate command errors and requires relatively less computational power among the others. However, performing calculations for the necessary derivatives of b-splines is a computationally expensive duty within the servo system sampling interval. Later, Cheng et. al.[12] focused on this problem and proposed a more efficient b-spline evaluation structure for faster computation. The structure reduced the total number of computations for b-spline derivatives by half. In this work, least squares approximation method is utilized to compress the G code blocks and to generate a smooth b-spline toolpath for five-axis machining. The toolpath is defined inverse kinematics independent enabling interpolation of the tool axis orientation vectors. Moreover, the generated tool path is smoothened, and its jerk content is reduced using linear and nonlinear optimization techniques to enhance the tracking performance of the servo system.  Chapter 2. Literature Review  12  2.4. Trajectory Generation for 3-Axis Motions Trajectory generation is vital for manufactured part quality. The generated trajectories must have continuous derivatives at least up to the second order to enable better tracking performance of the feed drive servo loop. Furthermore, acceleration and jerk profiles should have pre-specified limits in order to avoid saturating the actuators, or exciting unwanted dynamics in the servo or the mechanical structure. Various trajectory generation methods for 3 - Axis machining are summarized in the following. Prithschow [59] brought up jerk limited trajectories with trapezoidal or sine square acceleration profiles for better tracking. Makino and Ohde [53] used a universal Cam curve for jerk limited trajectory. Later, Simon and Isik [68] also used trigonometric functions for the same purpose. Similarly, Macfarlane and Croft [52] designed a jerk limited trajectory generation method where they utilized a sine wave approximation to determine the velocity ramp from zero acceleration to non-zero acceleration for robots. Especially, for 3 - Axis machining, quintic spline [78, 79, 80] and b-spline trajectory generation [38] methods are used, which enable acceleration and jerk continuity. Altintas and Erkorkmaz [18] used a similar trapezoidal acceleration profile to generate online jerk limited trajectories on quintic spline toolpaths. Particularly, jerk limited trajectories have high smoothness making it easier to track but also avoid high frequency harmonics in the acceleration profile that would result in degrading the tracking and contouring accuracy. Same researchers also used second quintic spline in their work to further smoothen the reference position commands interpolated at varying time period and reconstructed it at servo loop frequency. Pritschow and Rogers [60] implemented the resampling as micro decentralized interpolation within the advanced digital feed drives. Later, Altintas and Erkorkmaz [20] developed a trajectory generation method with cubic acceleration profiles for high speed milling operations. In this thesis, the fundamental requirements of the smooth trajectory such as acceleration and jerk continuity conditions are deployed to generate feedrate profiles for 5 - Axis machining. The research on generation trajectories for 5-Axis machine machining is presented with its critical issues in the following section.  Chapter 2. Literature Review  13  2.5. Critical Issues on 5-Axis Spline Trajectory Generation  The implementation of such C and C jerk limited spline trajectory techniques into 52  3  Axis is vital to adapt the process for HSM. Unlike to 3 - Axis, 5 - Axis machining centers have a non - cartesian kinematic structure [7,'64] that requires transformation of the tool tip position and tool axis orientations extracted by C A M systems into reference joint position commands using the inverse kinematics [51]. Some research has been done in robotic area on spline fitting and smooth trajectory generation techniques applied in the joint space [9, 39, 76]. Besides the robotic area, Chou and Yang [13, 14] presented an analytical study on the smooth command generation for 5 - Axis CNC machining by interpolating the joint angles with separate n order polynomials that are related to a single tool path parameter. In their th  work, they scheduled the curve parameter for smooth tool tip acceleration and constant velocity. Nevertheless, the kinematic limitations of machine tool's axes are not included in the study and can be violated because of tool path curvature and inverse kinematics transformation. Later, same researchers addressed that problem by using composite Hermite curves to generate piecewise constant velocity connected with smooth acceleration segments [15]. They included the machine tool's dynamic torque equation and axis jerk limitations to obtain the boundary conditions at the segment junctions and the maximum speed at the constant feed sections of the trajectory. Yang and Golub [85] also included a cutting forces predictor into the trajectory generation scheme and used inverse dynamics to enhance the tracking performance. On the other hand, these methods require a machine dependent post processing step at C A M systems to generate offline joint splines that must be updated for each kinematic structure. A more general and efficient method is pursued to interpolate the toolpath position and tool orientation (See Fig. 2.5) regardless of the machine tool and surface type. The primary obstacle is the interpolation of the orientation vectors, while ensuring that the entire interpolated orientation spline remains on the unit sphere. The standard mathematical tools to interpolate the tool orientation are based on homogeneous transformation matrices, Euler angles [51], Frenet frames [77, 82], ruled [37] and implicit surfaces [83]  as well as  quaternions. For instance, Siraivasan and Ge [71, 72] and Jutter [36] used quaternions to generate a dual three-space called the Image Space. They reflected the toolpath vectors into the image space and interpolated using rational Bezier and b-spline curves at C continuity. 3  Chapter 2. Literature Review  14  Cebula et. A l . [10] utilized the same idea using Nurbs and modified the orientation vectors in order to eliminate any singular location on the toolpath. / o = [o o o ]T x  y  2  Tool Orientation \ Vector  Figure 2.5: Tool Tip and Tool Orientation Vectors To tackle the problem of tool orientation interpolation, Lin and Koren [40, 41] implemented the ruled surface technique on 5 - Axis machining. They have first extracted a path curve on the part surface and used the normal or the ruling of the surface directly as the tool orientation vector during point and flank milling operations along the curve. Since the surface is already expressed in a parametric form, the tool orientation can be computed continuously for each parameter in real time by the NC unit. The tool tip velocity is then held constant by the aid of first order Taylor expansion. Spence et.al. [29] proposed a novel technique for smooth tool position and orientation interpolation using two nearly arc length parameterized quintic splines. The first quintic polynomial spline is used to interpolate the tool tip position vectors in cartesian coordinates. In addition to the tool tip spline, he defined a C quintic Bezier curve on the spherical 2  coordinates to be able to interpolate the unit tool axis orientation vectors (See Fig. 2.6). The spherical Bezier is also nearly arc length parameterized so that the constant tool rotation is achieved at equal curve parameter increment. He connected the travel distance to the tool orientation angle by the aid of a monotonic rational spline to achieve synchronized motion. This approach allowed generation of inverse kinematic independent trajectories for 5-Axis  Chapter 2. Literature Review  15  machining where the NC unit only applies transformation to compute the joint angles during real time interpolation.  a) Figure 2.6:  b) a) Position Spline — b) Orientation Spline  Due et. al. [47] addressed the tool orientation interpolation problem by using double bsplines defined practically in the workpiece coordinate system at the C A M stage. They presented the tool tip position vectors using a cubic b-spline. In addition, a second b-spline is employed in order to interpolate the position vectors of another point on the tool axis where the vector between the curves enables interpolation of continuous orientation vectors throughout the tool path. This method eliminates the need for a third curve definition for the coordination since both b-splines. are .parameterized by the same knot vector. Furthermore, the orientation of the tool can be modified using second b-spline to avoid singularities [1] or optimized [48] for flank milling operations. The double b-spline interpolation technique is also recognized by the Siemens NC units [69]. Parallel to the implementation of various spline interpolation techniques, development of feed optimal trajectory is necessary for five-axis high speed machining have gained high importance. It is a challenge to generate a smooth trajectory that enables the spline tool path to be traveled in minimum time without violating the physical limits of the axes drives. A significant amount of research has been done on time optimal trajectory generation for robots [6] and Cartesian machine tools [8]. However, those trajectories considered the maximum  16  Chapter 2. Literature Review  torque or accelerations limitations resulting in a non-smooth trajectory that leads higher contouring errors. Constantinescu and Croft [16] used joint variables expressed as cubic polynomials for continuity. They optimized the cubic spline trajectory with respect to the torque rates for minimum cycle time. Results have showed enhanced tracking performance. On the other hand, the implementation of such optimization architecture is not adequate for machine tools' NC units because of its computational expense. Week et. al.[81] addressed that problem and developed a practical method to be used for machine tools processing spline toolpaths. He considered the worst curvature case in each spline segment to determine the maximum feedrate for that segment and connected the piecewise constant feed segments together with smooth accelerations transients. Although his method was easy to implement, higher velocities would be achieved in non-critical sections of the segments by using a continuously varying smooth feed profiles instead of piecewise constant. If continuously varying feed profiles are preferred to use instead of piecewise constant velocity segments, the difficult task is obtaining a smooth velocity profile, which can tracked with less error by the drives. Flash and Hogan [28] first suggested minimum jerk trajectories for describing the smoothness of unconstrained human arm movements. Later, Kyriakopoulis and Saridis [44] proved that tracking performance is enhanced by minimizing the jerk content and implemented minimum jerk trajectories for robotic manipulators. Later they, focused on practically optimizing the trajectory for both minimum jerk and cycle time and minimized the maximum jerk value of the trajectory [45]. Piazzi and Visioli [57] also worked on minimizing the maximum magnitude of jerk using a global optimization algorithm for cubic trajectories by modulating the time durations of spline segments. In addition, Simon and Isik [68] used trigonometric splines for minimum jerk trajectory. They showed that the problem can be reduced into quadratic minimization problem, which makes it possible for practical implementation.  ,  Altintas and Erkorkmaz [3] developed an efficient trajectory generation method for quintic spline tool paths, where the vitality of minimum jerk feed profile is combined with the cycle time objective. In order to schedule the feedrate continuously on the toolpath, they defined the displacement profile as a C quintic spline of time and optimized its coefficients 3  for minimum jerk. The minimum cycle time objective is then injected by modulating the time  Chapter 2. Literature Review  17  durations of each spline segment while respecting the machine tool's torque, velocity and jerk limits using an iterative algorithm. The algorithm is applied with a windowing technique for efficient processing and real time application. A b-spline toolpath generation algorithm has developed in this thesis for 5-Axis machining based on the fitting a double spline toolpath on discrete data. The position and tool orientation errors are limited by user specified tolerances. Furthermore, a minimum jerk, minimum time trajectory is generated optimizing the quintic spline feed profile with respect to the machine tool's axis drive limits for the b-spline toolpath. The generated trajectory is then interpolated using the feed correction polynomial idea in real time.  75  Chapter 3 Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining 3.1. Introduction  Reference trajectory generation plays a key role in the computer control of machine tools. Generated trajectories must have smooth'kinematic profiles in order to maintain high tracking accuracy, and to avoid exciting the natural modes of the mechanical structure or servo control system. Especially for five-axis machining, the generation of a smooth trajectory is vital for the part quality. For five-axis machining, after the design of the part, the primary machining toolpath is sampled for linear or circular interpolation by C A M software. For instance, these sampling algorithms extract relatively large number of points from the toolpath where the rate of change of curvature change is high and fewer points where the rate of change of curvature is low. This results in a G code containing thousands of polygonal segments to machine a sculptured surface. Feeding the machine tool with large number of discrete points leads to discontinuous movements of the axes for high speed machining since the line segmentation does not ensure continuity of the kinematic profiles. Therefore, spline trajectory generation techniques have become widely adopted in machining aerospace parts, dies, and molds. Their main advantage is that they provide more continuous feed motion than multiple linear or circular segments, and result in shorter machining time as well as better surface geometry. The objective of the work presented in this chapter is to generate spline toolpaths for five-axis machine tools from initially sampled discrete toolpath data. The desired spline toolpaths are defined in terms of double b-splines [47]. The position and orientation of the tool can be interpolated directly using the toolpath format. Furthermore, two optimization techniques are presented in order to smooth the motion of a five-axis machine tool following the b-spline toolpath. Herein, a b-spline toolpath fitting algorithm is developed by utilizing constrained least squares techniques in the toolpath generation part. It accepts discrete toolpath data consisting  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  19  of Tool Center Points (TCP) with tool orientations and generates two b-splines for a specified machining tolerance. The first b-spline is fitted to the translational component of the toolpath. It is used for the interpolation of the tool tip position. The second B-Spline fitted to the second set of points is used to interpolate the orientation component of the trajectory. In the following optimization work, the quality of the b-spline tool path is improved. For smoother motion of the tool tip, the third derivative (i.e. Jerk) minimization objective is inserted into the least squares algorithm that delivers a smoother geometric profile of the tool tip to the machine tool. In addition, a nonlinear optimization algorithm is proposed for the upper b-spline path. This algorithm modifies the upper b-spline path by penalizing the angular jerk component of the toolpath. The resultant geometrical toolpath will enhance the tracking performance of the machine tool.  CAD Model I  Discretization Tool Orientation Data Tool Center Points —  —k Generate Tool Tip (Lower) B-spline  i i1 *>  Optimize the Tool Tip spline  1  y t Generate Second (Upper) B-spline for the Tool Orientation  j [ i  Generation of the J Double B-spline Tooj-pathj  Optimize the Tool Orientation spline  i  [ j  Minumum Jerk | spline Tool-path !  Figure 3.1: Flow Chart of Optimal B-spline Toolpath Generation Firstly, a brief background is given for the b-spline representations and computations followed by the double b-spline polynomial representation of the toolpath in Section 3.2. The least squares based fitting of the double b-spline toolpath is presented in Section 3.3. The optimization algorithms applied on the spline toolpath are presented in the last section.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  20  3.2. B-Spline Background  Cubic and quintic splines are being widely used for NC toolpath generation [78, 79, 80] as a replacement of line segmentation. For smoother trajectories, cubic polynomial splines are employed since they can ensure C continuity allowing interpolation of acceleration 2  continuous kinematic profiles. Similarly, quintic splines are utilized in the toolpath generation for a more consistent parameterization [78]. Unlike cubic and quintic polynomial splines, b-splines have only recently been adapted into the NC toolpath generation for machining because of their higher computational cost. However, b-splines overcome various deficiencies of previous splining techniques as follows:  •  B-splines are compatible with the CAD Systems, since most parts are designed using surfaces generated by b-splines.  •  B-splines are smoother curves.  •  It is easier to modify and predict the shape of the curve with the aid of control points.  •  The curve can be modified locally without affecting the shape of the complete spline.  A b-spline is illustrated in Fig. 3.2. b-splines are a generalization of Bezier curves. They are defined by a control polygon consisting of n +1 control points P  n  the shape of the curve.  i = 0,... n that form  21  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  y Control Points  p+1  Knot Vector  1  ^+1  '  Figure 3.2: B-spline Toolpath By  definition, a p  th  degree b-spline constructed between points P  0  and P  n  is  parameterized in the interval of [0 l] by a knot vector, U. The knot vector contains m = n + p + \ knots (//) that are initially allocated spline parameters (w) for certain positions on the b-spline (See Fig. 3.2). The knot vector of the b-spline is given as:  Repetitive p+1 Knots  Repetitive p+1 Knots  where  (3.1)  Mj < fi ,j j+i  = 0,1,...,m, ps[0  The consecutive knots [jUj p  J+l  1].  ] in the knot vector are called as knot spans. In addition,  the knot vector contains p + l repeated knots at the beginning and end of the knot vector in  Chapter 3. Spline Toolpath GenerationfromDiscrete Toolpath Data for 5-Axis Machining order to anchor the b-spline at those control points. For instance, the knot vector of a cubic (p=3) b-spline is given as:  U = [0,0,0,0,p ,p ,...pj,...,//„_!,fi„,U,U ] • p+i  p+2  The evaluation of b-spline is a computationally expensive. Any position on cubic or quintic splines are expressed by a single explicit polynomial equation such as: s (u) = a y  + a _y-  +...+a u  [  x  x  0  xn  xS)  S (u) = a u" + a _ u ~ +... + a u°. n  y  yn  yn  x  x  yQ  However, since b-splines are built up by intermittent piecewise polynomials (basis functions) that blend the control points. A recursive algorithm is used in order to compute the basis functions. The i b-spline basis function of degree p is denoted by N . For example th  ip  a zero degree basis functions is defined as:  ( )=\ ^- ~^M Uf  N  U  1,0  U  ( 3 2 )  0 otherwise.  As indicated in Eq.(3.2), zero degree basis functions are step functions. Higher degree basis functions [N ,N ,...,N ] iX  i2  ip  are linear combinations of two lower degree basis  functions expressed according to the following recursive scheme:  22  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining M,  ~u  +P+l  MM  ~Mi  Mi+p+\ — MM  M«)=^-",>)+ MM  23  -M>  M \  ~u  i+P+  Mi+p+\ — MM  (3.3) M x  -u  i+P+  MM  Ni. (u) P  =  -Mi  Mi+ +i  M  P  ~u  i+P+i  ^ ^ N ^ ( u ) MM  NM. M  ~MM  P  +  Mi+ +i  -M  P  MM  and various basis functions are illustrated in Fig. 3.3. a) Step Basis Functions ^4,0  *3,0  [u] ~  r  ~3  ~  ~3  Ti  ~  ih  "  b) 1st Degree Basis Functions  [u] H2  c) 2  4  t!  0  1  3  ^2  H  1  0  3  H2  H  H2  H3  3  Degree Basis Functions  nd  N , 0  2  I  >.j 0  H  1  Hi  M2  1  !-».*  H3  1  0  L_ Hi  >>4 H2  I  H3  U*  1  0  Hi  1  Figure 3.3: B-spline Basis Functions  Once the non-zero basis functions are computed for the spline parameter (w), a position on the b-spline is obtained by the sum of only the p degree basis functions and control th  points as:  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining P(u) = N^ {u)-P +N {u)-P p  0  hp  24  + ... + N _ ^  v  n  x  (3.4)  Similarly, derivatives of the b-spline curves are calculated using the derivatives of basis functions. The first derivative of the p  th  degree basis function with respect to the curve  parameter is given by a similar recursive scheme:  dw  (3.5)  — ^ w » .  M -Mi  Mi+p+i  l+P  Mi+i  and the k derivative of the basis functions is expressed as: th  ) > =P  du (k)  Mi+p+i  (3.6)  Mi+i  Then, the derivative of the b-spline is obtained by:  d<M")_ X>) „ d(  d^N„_ (u)  du  du  lp  du  ^d x>)  d^N (u)  (k  hp  d^N , {u) n p  du  (3.7)  p  du (k)  3.3. Double B-Spline Polynomial Representation A toolpath consists of the tool tip position and associated orientations of the tool axis in five-axis machining. In order to maintain this information of the tool motion, a spline toolpath format is put forward by Lartigue at. al.[47]. The proposed toolpath consists of two splines, which are dedicated for the following:  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  •  25  The locus of the tool tip position (x',y',z') is given by a spline, P(u) in part coordinates.  •  Another point on the tool axis (x ,y ,z ) h  h  h  is represented by a second spline,  Q{u).  Figure 3.4: Double B-Spline Toolpath Representation The overall representation of the spline toolpath is illustrated in Fig. 3.4. The lower spline P(u) is used to express only the tool tip position as in regular 3-Axis machining and the upper spline Q{u) interpolates another point on the tool axis. The relative distance between P(u) and Q{u) is defined as H, which can be defined depending on the machining operation such as flank or point milling. For flank milling, the relative distance refers to the metal cutting contact height that corresponds to the active part of the tool or depth of cut. Similarly, for point milling H is the diameter of the tool's ballend as illustrated in Fig. 3.5.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  26  Tool Axis(A(u))  Control Points(P ) h  I  Control Points(P ) t  Figure 3.5: Upper and Lower B-splines The splines of the toolpath, P(u) and Q{u) are given in non-uniform b-spline form. The expression of lower and upper b-splines defined on the same knot vector, (U) is expressed as:  x'{u) P(u) = y'(u)  u eU /=<)  1=0  z'(u)_  (3.8)  x (u) h  fix  Q(u) = y (u) h  z {u)  where  P  T  =[P  X  P  Y  Qy  1=0  ,<=0  h  PJ Z  and Q =[Q t  X  Q  y  Q f, 2  ueU  i = 0,..,n  are the n + 1 set control  points of the lower and upper b-splines, respectively. The degree of the b-splines is selected either p=3 or higher depending on the continuity requirement. For instance, the cubic double B-Spline toolpath allows definition of a (c ) continuous motion. Similarly, a 5 degree 2  th  (quintic) b-spline toolpath enables interpolation of jerk continuous (c ) motion. 4  The tool tip position is directly interpolated by P(u). However, in order to interpolate the tool orientation, both b-splines are utilized. As shown in Fig. 3.5, the vector A(u) defines the  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  27  tool axis. It can be expressed for each spline parameter simply by subtracting the lower bspline from the upper one as:  ~x -x'~ h  X ~  =  A(u) = Q{u)-P{u) = y -y' h  z -z'  (3.9)  = A i + A j + A k.  >  x  A  h  y  t  A .  The unit orientation vector of the tool axis,0(w) is computed similarly when A(u) is normalized by the relative distance between P{u) and Q(u) as:  r>'i + O j +r,7 O O k =x(")' ^ ^i A  0(«) =  x  v  z  =  y( )-+ ^ - ^ j 4+f X ^ > 7^ k  A  u  (3.10)  r  o.  where 0(u) is also called orientation spline of the trajectory. It is expressed in terms of the upper and lower b-spline control points by:  ix»-fe - i) p  0{u)  = f^  ir,Mec7.  (3.11)  1=0  Hence, the objective of continuously expressing the tool axis orientation and tool tip position is achieved during simultaneous interpolation of both splines using Eq. (3.8).  3.4. Generation of Lower B-Spline  The tool tip positions are initially given as a batch of discrete tool center points (TCP) by C A M systems. The fitting algorithm for generation of the lower b-spline through TCP data is presented in this section. For 5 - Axis machining of free-form surfaces, large number of TCPs is needed and interpolating all the data may result in a non-smooth spline. In order to  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  28  generate a smooth b-spline toolpath, a least squares approximation is used instead of the exact interpolation. A tolerance band is defined to limit the least squares fitting errors to respect the machining tolerances. Furthermore, the algorithm is extended as a general least squares to impose any specific position, tangent and higher derivative constraints onto the toolpath. The lower b-spline is already defined in Eq.(3.8) as,  x'{u) u e U,  P(u) = z'(u)  and the above expression shows that modification of a control point affects only the associated axes components of the b-spline. Therefore, the least squares can be applied to each axis component separately. The formulation used for fitting the x' component of the lower b-spline is presented here and same formulation is then applied to the rest of the bspline components.  B-spline  u=u =o 0  Figure  3.6: Initial Spline Parameter Calculation  The tool tip position (q' ) is given by M + l data points, q' q[ k  q[ [q'x =  k  q[, q'  M  where  Vy, <lz J • The least squares algorithm starts with seeking the best prediction, k  q' for any tool tip point q' at x'(u\ x  0>  x  _ on the lower b-spline,  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  <7 = fc„ k\  C  where u = [w  0  2  u  {  -  q'J  u  ...  2  ^x {uhN {u)-P N Ju)-P +... l  QtP  u]  x0+  l  xl  +  29  N Ju)-P „(^2) n  x  are the parameter values on the b-spline corresponding to  M  the given tool center points (q' ) shown in Fig. 3.6. The chord length parameterization k  method [58] is applied to calculate the initial spline parameters. This method is used for most of the data distributions and gives a rough estimation of the spline parameters while only considering the Euclidian distance between the points. The total cord length between the given data points is defined as,  S = t\q' c  k  k=\  ~ql\  =t M - < J  +  ( < ~ < J  +  ( < - < J  •  <" ) 3  13  k=\  The initial and final parameter values are set as u = 0 , u = 1 and the internal spline 0  M  parameters are allocated proportional to the consecutive chord lengths as:  u =u _ +k  k  k = l,...M-l.  t  (3.14)  Thereafter, the knot vector, (U) of the b-spline is required for the parameterization and setting up the equations for the least squares fitting. The initial and final repetitive knots in the knot vector are defined by Eq.(3.1) initially as:  Mo= — = Mo=°  a  n  d  M„  +l  =--- = / W i  =  1  -  The remaining knots are obtained as described by Piegl [58]. Defining d as a positive real number, denoted by / = int(^) the largest integer such that  i<d,  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining m+l n-p + l  30  (3.15)  where n is the number of control points, p is the degree of the b-spline and m = n + p +1 The internal knots in the knot vector are calculated by:  i = int(j-d) (3.16)  a - j -d — i M  = (l -  P+j  j =  h...,n-p  Using Eq.(3.16) and Eq.(3.1) the knot vector of the spline toolpath is formed. Corresponding position predictions for the data points are obtained on the b-spline using Eq. (3.14), which results in the following system of equations:  ^I,P("O)  *;  =  •  ^0,„(« ) 2  N (u ) np  (3.17)  2  _^0.p("«)  where 0 is a M +1 x « +1 matrix consist of basis functions (iV, («)) evaluated at each u , k  and P is the unknown control point vector. The errors (residues) of the predictions are x  defined as the deviation of the approximated points q' = x' (w)| _ on the b-spline from the Xt  actual tool tip positions, -q' expressed as:  e -q[-q' =q' -Y^N, {u)P =q' -0P x  x  x  p  xi  x  X  /=o  The objective function of the least squares to be minimized is then;  (3.18)  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  1 1 2 " " 2 T  <?: -t J") >] N  P  U-£Nju)P }  = Uq' -0P  x  x  31  ) (q' -0P T  X  x  X  ). (3.19)  The formulation of least squares implies that the derivative of the objective function with respect to the control points P is set to zero. x  dJ,  =  dP  -0 (q' -0P ) T  x  0 0P  (3.20)  =0q[  T  X  P  x  ^0] 0q' X  x  x  Eq.(3.20) gives the control points for a regular least squares fit. The b-spline can also be constrained at some certain points. For instance, position constraints at tool tip points can be introduced in order to eliminate any fitting errors. In addition, tangents or possibly higher derivatives required to be imposed on the curve. In order to enforce the desired shape constraints, a general linear system of equations is composed in the constrained least squares algorithm. For the specific position constraints x'(w)| _ = 4  (o)  ,k = 0...l , linear equations are 0  grouped as follows: N , ( o) u  0 P  i ( o)  N  N {u )  u  P  2p  0  ...  N„ {u ) p  0  (3.21) (o)  3  N , k) 0 P  tf.,(«J  2, k)  N  P  -  ». k)  N  P  L<°>  Tangent,  dx'(w) du  d x'{u) 2  ' =$\  k = 0...l  1  and  normal  du  2  =£ \  k = 0...l  2  constraints can be further imposed in the following matrix equations (Eq.(3.7)) by using the  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  derivatives  of  the  dN  basis  functions  given  as  N (u ) jpu  =  k  32  AN;  and  I.P  du  2  du  I.P  , which are calculated by Eq. (3.6),  u M  N {u )  2  ho (1)  N  Q  ,  P  XO.PM)  lpu  N {u )  ...  2, M)  -  0  2pu  ^.PM)  0  n  p  N {u ) npu  N  0  »,A)  (3.22) 0 .  N  4,  P  u u M  0,p«u(u, )  N  2  l p » u M  N  l uu{u, )  2.puu(u, )  N  P  •••  2.puuM  N  N  2  •••  2  N  N  n . p u u M  n . p  U  U  ( » ,  2  P.  )  (oy  L  The resultant constraint matrix is L =  where r = (/„ +1)+ (/, +1)+ (l +1) and the 2  Jrxn+1  overall constraint vector containing position and derivative constraints is £ = £(2)  Hence, the lower b-spline is obtained by solving the following optimization problem:  J = m i n i f e - 0P J(q' - 0P ) L  X  x  X  Subject to:  (3.23)  L-P =£, x  which is a linear quadratic minimization problem, and can be solved using Lagrange Multipliers method. Introducing the Lagrange multipliers: A - [4> augmented objective function is written as:  4  Ai  '" 4-]>  t n e  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  J , {P,A) = i f e - 0P J(q' L  X  x  Equating the partial derivatives to zero  - 0P )+ A(L  x  yields the linear equation  dA  dP  (3.24)  -P -$)  X  33  J  system:  0 0P T  LP,  +L A = 0 q' T  X  T  ~0 0 T  x  L  =£  L~  ~0 q' ~  T  0  T  x  A  From Eq.(3.25), a solution can be obtained only if 0-0  T  singular. Every knot span, [Pj,p \ j+X  condition of 0-0 . T  (3.25)  is positive-definite and non-  must contain one u to satisfy the non-singularity  Eq.(3.16) guarantees that every knot span contains at least one u ,  which eliminates the possibility of singularity. Furthermore, provided that the constraints are linearly independent and less than the unknowns (r < n), a unique solution is obtained for the constrained least squares fit.  3.4.1. Tolerance Definition and Parameter Correction  The general least square algorithms are simply based on minimizing the total vertical errors in each axis component of the b-spline. However, in this application, it must be ensured that the approximated b-spline stays within a specified accuracy, namely the machining tolerance. A fitting error limit (e ) is defined to control deviation (e) on the max  lower b-spline from the reference tool center points. In Fig. 3.7, the tool center points, the proposed b-spline fit and the definition of the tolerance band to limit the errors can be seen. The tolerance band is defined as a pipe with a radius of (e ) in 3 dimensional space max  covering the fitted b-spline.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  34  6 Error (e) Figure 3.7: Tolerance Definition for Lower B-spline Generation  The deviation on the lower b-spline from every tool center point is required to check the errors. Corresponding spline parameters («) to the tool center points [q') are already obtained by the chord length approximation. However, as shown in Fig. 3.8, this approximation does not ensure the correct spline parameter for each tool center point on the generated b-spline. The closest spline parameters must be identified in order to calculate accurate norm of deviation from given points, which is done by scanning throughout b-spline for the correct parameter. The vector between the lower b-spline P = \x §' - [q'x q'  y'  z'J  and the given tool tip point  q'z J * written as: s  y  D = (x<-q' ) x  +  (y>-q')j  +  (z'-q:yc.  (3.26)  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  35  It is known that the vector D must be perpendicular to tangent of the spline at u for the closest point. A dot product function is formed to ensure the condition as shown in Fig. 3.8; f( )=«M\ (u)- '). u  P  Figure 3.8:  (3.27)  q  Spline Parameter Correction  The root of the Eq. (3.27) is found by Newton Raphson iteration. Parameter w. is obtained by the i' iteration, h  dP(u) -u.  /(«)  • u, -  du d P{u)  du  (P(u)-a>)  2  du  2  (P(u)-q')  +  (3.28)  d/>(«)/ dw x  The convergence of the iteration in Eq. (3.28) to the correct spline parameter depends on the initial guess. A very close initial guess is simply thew values computed initially by the chord length approximation by Eq. (3.14). Two zero tolerances, shown on Fig. 3.9, are defined to indicate the convergence: s, is the measure of the Euclidean distance and s is to indicate zero cosine value. 2  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  a)  jq  Measure of Zero Cosinefed  t  36  Tool - Tip Point  D  Accurate Spline Parameter b)  Measure of Euclidean Distance^)  Tool - Tip Point Accurate Spline Parameter Figure 3.9: Zero Tolerance Definitions a) Zero Cosine — b) Point Coincidence  Using the above zero tolerances, the criteria to terminate the iteration is as follows:  d)if  (3.29)  P(4 = -<l' u  u  the point lies on the b-spline at the calculated parameter (point coincidence) or,  (2)  if  dP{u) (P(u)-q') du dP{u) (Pi")-*'] du  (3.30)  the vector D is the normal of the curve at the spline parameter and condition for Eq.(3.27) is satisfied or,  (3)  if  dP{u) du .  (3.31)  37  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  the spline parameter does not change significantly during the Newton Raphson iteration (Eq. (3.28)). For basic machining operations, s  x  in Eq. (3.29) and Eq. (3.31) can be taken as  1(T [mm] to indicate that the point lies on the b-spline at the searched spline parameter. On 6  the other hand, for precision machining this threshold value should be selected more conservatively. The zero cosine measure, e in Eq.(3.30) is selected less than e such as 2  x  1(T . 8  If any of the above conditions is satisfied (Eq. (3.29), Eq. (3.30) or Eq. (3.31) ), the iteration is terminated and correct norm of the error is found at u, as:  Dl  Increasing the number of control points, n improves the control on the spline shape. Eventually, the maximum norm of the error can be reduced to the desired value by increasing the total number of control points in an iterative approach during least squares fit. The procedure is shown in Fig. 3.10.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  38  "F  ^3 - Tool Center Points (q ) l  0-  , - Max Position Tolerance (e )| max  Control Points n = M+1 / 20 Arc Length Parameterization for Initial Knot Vector (U) Insert an extra Knotf-H  n  Find the (knot span (i) of ( e 1—(n = n + l)-<-  Fit Lower B-spline P(u) 1  max  pheck errors (e)J  j No  if e<e„  Yes  End  Figure 3.10: Lower B-spline Fitting Algorithm  The fitting process of the lower b-spline starts with less number of control points. The number of initial control points can be judged for free-form surface machining as 5[%] of the total given TCP (n = ( M + l)/20). Using Eq. (3.1), Eq. (3.16), the knot vector is generated and the lower b-spline is fitted using Eq. (3.25). Following that, resultant errors are checked using Eq. (3.27) and Eq. (3.32). The knot span, where the maximum tolerance violation occurs, is recorded as un-converged. An extra knot, ju is added in the middle of the unconverged knot span by averaging the neighbouring knots as:  Mj+M  J+i  (3.33)  Adding extra knots only into the specific knot spans will assist reducing the errors focused mainly at this portion of the b-spline. The b-spline is then fitted again with the new knot vector using Eq. (3.25). Errors are checked after each fit until all the knot spans are converged and the max error is kept under given limit.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  39  3.5. Generation of the Upper Spline  The fitting of P{u) through the set of tool center points in the previous section. Similarly, a second set of points on the tool axis is required to generate the upper b-spline Q(u). Firstly, the reference tool center points are updated on the lower b-spline using spline parameters obtained from Eq. (3.28). The updated points are offset in the direction of reference orientation vectors O = [O r  O  r  O ] by the amount of H as,  r  x  r  q =P{u )+H-O h  h  k  (3.34)  k = 0,...,M  r  k  where q  T  z  k  is the new set of reference points in cartesian coordinates for another location on  the tool axis. The upper b-spline Q(u) is then fitted through q using Eq. (3.25) while h  checking only for the orientation errors e . The orientation error is defined as angle between e  given reference orientation O and the approximated orientation, O as: r  . /  \  e =tan"  \O xO \ r  k  f  s  i  e  k  n  ( J  e  k  where  o -O r  k  t  e  ( e ),  cos  ,  ,  - tan"  (3.35)  <Z-o  K  kj  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  40  An angular tolerance cone is illustrated in Fig. 3.11 to indicate the limits of tool orientation errors. Similar to the upper b-spline generation, an iterative fitting process is applied until the maximum tool orientation error is below e  e  1 2 3 4 5 6  .  Upper B-Spline Data Points (q ) Upper B-Spline (Q(u)) Angular Tolerance Cone Orientation Error(e0) Tool Orientation(O) Reference Orientation(O ) h  r  Figure 3.11: Upper B-spline Generation and Tolerance Cone Definition for Orientation  The algorithm shown in Fig. 3.12 summarizes the least squares fitting of the double bspline toolpath. The resulting knot vector of the upper spline might have more knots then the lower one or the vice versa. To make the trajectories synchronized, both of the knot vectors must be equal. Therefore, abundant knots are inserted into the knot vector of the b-splines. The method by Piegl [58] can be used to insert sparse knots into the knot vectors. Added redundant control points do not change the shape of the curve or the parameterization, and the resultant knot vector will preserve desired synchronization.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  | - Tool Center Points (q) 3 - Orientation Vectors (O ) n - M a x Position Tolerance (e, max/ ~7- M a x Angular Tolerance (eg, max. — - B-Spline Degree (p) r  Control Points n = M+1 / 2 0 A r c Length Parameterization for Initial Knot Vector (U)  r  Insert an /" extra Knot—W  Find the knot s p a n (i) of ( e t_(n = n  +7)4-  Fit Lower B-spline P(u) C h e c k position errors (e)  m a x  j No  if e < e  m a x  1 Yes U s e s a m e knot vector (U)  Q(u) C h e c k rotation errors (eg)  'f q<e0n  >  Insert a n lextra Knot (m)  Fit UpperB-spline  No  1 Yes  | _  o ^  3Xn+1  3Xn+1 -U  n  >(n = n + 1 ) J  Insert abundant knots to generate single knot vector  b - P l Q_ -Qi  Q  Find the knot s p a n (i) of (e0  1Xn+p+l|  Figure 3.12: Double B-spline Generation Algorithm  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  42  3.6. Smoothing the Double Spline Toolpath  3.6.1. Introduction  The quality of the fitted curves is crucial for high speed machining. When machining at very high feedrates, the geometrical properties of the curved toolpath directly affects the tracking performance of machine tool. For an improved tracking of the toolpath, the resultant trajectory should deliver very smooth kinematic profiles to the drives. The dominant factor affecting the tracking performance is identified as the jerk component [28] of the trajectory. The total jerk depends on both the geometric properties of the toolpath such as the change of the curvature, as well as the applied kinematic profiles that build up the complete trajectory. For instance, if the geometric profile shows sharp changes and high roughness, a varying or very conservative feedrate must be applied in the trajectory generation to avoid contouring errors and possible kinematic limit violations of the drives causing possible actuator saturations. Therefore, the geometric toolpath is desired to be as smooth as possible containing minimum changes in the curvature i.e. jerk components. This allows imposing more constant and smooth feed sections in the kinematic profile at the trajectory generation step. As a result, the total machining time will be reduced and the tracking performance of the drives is improved. Geometric jerk is the third derivative of the position vector, P(u) with respect to the arc length (s). Minimization of a cost function based on geometric jerk can be obtained through nonlinear  optimization methods.  However, this process  introduces  a significant  computational load to the trajectory generation. In [22] it was proposed that minimizing the third derivative of the position vector with respect to the spline parameter (w) results in significantly reduced jerk profiles as well. Therefore, this approach is adapted in order to minimize the jerk content of the trajectory. As mentioned, the motivation of this section is to generate a smoother geometric profile of the toolpath, which is easier to track. Two optimization algorithms are presented to modify the double b-spline toolpath. Firstly, the lower b-spline is addressed. Secondly, the orientation of the tool is optimized by modifying the upper spline for smoother rotation of the tool.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  43  3.6.2. Smoothing the Lower B-Spline  In the linear least squares fitting process, the quality of the curves mainly depends on the constructed knot vector and the control point positions. The initial knot vector is created by the chord length method that might cause oscillations or twists on the b-spline [56]. Moreover, during the least squares fitting process of the b-spline toolpath, extra knots are added into the knot vector (see Fig. 3.12). This can distort the knot distribution and may cause severe fluctuations on the b-spline derivatives. These shortcomings of the least squares approximation are addressed, and the overall smoothness of the toolpath is improved in the following using the jerk minimization objective. The objective in this section is to smooth the lower b-spline of the toolpath by a linear optimization algorithm, which mainly contributes to the motion of cartesian axes of the machine tool. The tool tip position is given by the lower b-spline as P(u) = [x' (u) y' («)  z (u)f . The  geometric jerk is defined as the third derivative of P(u) with respect to the spline parameter and its minimization can be posed as an integrated square minimization problem as;  J =mm[(pl-P y . K  mu  (3.36)  u  d, is the vector of b-spline control points, ( \P  X  the lower b-spline calculates as,  P  Y  P^)  and P  UUU  is the third derivative of  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  P.  = y  where  uuu  Ni+lp-3to  N  {M,  P  =  dw  P(P~ XP- ) 1  to  l+2.p- («)  N  3  ~ MM  i+p2  i+P  [Pi* -MilMi+p-t -Mi) J  Z+l.p-3  M  M -2 ~ MM  Mnp-3 ~ Mi  N,  44  M -i ~ MM i+P  -MilMi  'Mi)  +P  + p  2  3  N Ju)  Mi+p-i-p,  MI+P-I-MM  jy (»)  N Ju)  i+hp  [MI+P+I — Mi+\\Mi+  p  N _ (u)  MtP  MtP  ;  +P  3  M -M  Mi -i-Mi 2  i+p  +  M  [MI+P+I — Mi+i^Mi+p+i — Mi+2)  — Mi+\)  (3.37)  Expressed in Eq.(3.37), the axis control points [[P P x  y  P \) affect only the third z  derivatives of associated axis components. For instance, modifying the control point P  x  would change the derivatives of only x axis components. Therefore, P  uuu  in each direction  can be minimized separately. The jerk objective function J in Eq. (3.36) is then expressed K  separately for each axis component as follows:  J k  ^ ^ " l  [( uuu-Puuu)du P  PPP« where x' , y' , y' um  um  um  the J , J , J K  K  K  <=>  '  =  and  m  P  i n  J  i ( "»» ^ " ' X  d  *>  J  =min|(z^„) d  =  m  P  i n  1 ^'"" ^ " ' u  d  (3.38)  2  Kz  M  are the axis derivative components (Eq.(3.37)) of lower b-spline and  are the objective functions to be minimized for each axis component,  respectively. The solution process of minimizing the jerk objective function J  K  for the JC  component is shown in the following. The integral square of the third derivative (jerk) throughout the b-spline (« e [0 l]) is written as,  J  K t  - min ^[x'uuJ  • x' )du . uuu  (3.39)  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  Substituting the expression for x'  from Eq. (3.37), the above integral is rewritten as:  uuu  JK, =(P ) (f  where  N  = [N  UUU  )=(P ) {KtP )  -K )dujp  r  X  m  (3.40)  T  X  x  N {u)  (U)  x  N {u)  lpum  0PUUU  45  ... N {u)]  and K is a n + \xn + l  NPuuu  Xpuuu  symmetric weight matrix, which is only a function of the integral of basis function derivatives. For any given knot distribution, K matrix is constructed for the b-spline parameter interval u e [0 l] as:  | o, N  P  uuu ( « K .  P  uuu (")d"  I ^o,p uuu ( " K p l N  Up  m  u  (u)N  LP  K  uuu  | N  „„ (")d" [uYu  | N  I  P  uuu  0  P  u  [u)N  m  2P  {u)N  2P  uuu  u  m  {u)du  (u)du  l 2. uuu(u)N , uuu(u)du N  P  2 P  SYM t 0.puuu{u)N„, {u)du N  puuu  t lpuMN„, (u)du N  pmu  i 2.pu„u{u)N , {u)Au N  • (3-41)  n pmu  l n,puM^ . {u)du N  n pmu  When the basis functions are in Bezier form, the weight matrix (K) can be expressed analytically since the knot vector is uniform and basis functions are repetitive. Identifying an analytical form for a single basis function will be enough to compute the integrals. On the other hand, during the fitting process of non-uniform b-splines, the knot vector changes at each iteration when extra knots are added. It is computationally expensive to calculate each basis function analytically during the least squares fitting. Therefore, the K matrix is computed by numerical integration methods such as using Simpson's 3/8 rule, which is relatively efficient.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  46  The previous optimization problem for the least squares fitting was given by Eq. (3.25) as  A  (') P A  =\ k -  ®* J  -  p  -P -Z),  )+ A(L  X  and it is modified by adding the minimum jerk objective function defined in Eq. (3.40) for the smoothing action. The augmented minimization problem can be represented as:  minfj, +J .l=min— r(q' -*P i  x  x  ffc-®P )+a(P X  x  Fitting Affect  }{KiP )  (3.42)  x  Smootheing Affect  Subject to: L • P = % x  which is also a linear quadratic minimization problem, and can be solved by Introducing the Lagrange multipliers: A = [AQ 4  4  4]> the resultant linear system of equation is  obtained as,  (0 0 T  r  + aK)P  + L A = y0 q' T  x  L~  T  'y0 q'  T  x  0  L-P =$ x  T  x  A  (3.43)  4  where a and y are the weights of the multi-objective minimization problem, a is the weight for smoothing effect and y specifies the fitting effect. Increasing y would increase the fitting affect. Therefore, the overall errors will decrease but the smoothness of the resultant spline would decline. Respectively, a increases the smoothing affect during fitting. The two parameters must be tuned for the fitting requirement. By experience, the parameters can be selected around^ = 10 ~10 and « = 10" ~10" . where their affect is investigated during 6  simulations in Section 3.7.  7  3  2  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  47  3.6.3. Angular Jerk Minimization  In a similar manner to the optimization objective mentioned in the previous section, the angular jerk profile throughout the toolpath can be augmented. Similarly, angular Jerk is defined as the third derivative of the orientation spline O (u) = d 0(w)/du . Minimizing the 3  3  uuu  3 derivative of the orientation spline helps smoothing the rotational motion of the tool. rd  Upper Control Points(P ) h  Figure 3.13: Tool Axis Vector The orientation spline was defined in section 3.3 as,  (3.44)  0(u)  where the tool axis vector A = A(u) shown in Fig. 3.13 is defined between the upper and lower b-splines of the toolpath as, (  '< A{u) =  >  A  A  z  = Q{u)-P(u)  1  - Qy  = J Nijuf j  i=0  \\  >," P,  JJ  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  where [P  P  X  P ] and [Q  Q  r  Y  Z  X  48  Q ] are the control points of lower and the upper br  Y  Z  splines, respectively. In order to obtain the derivatives of 0(u), the norm of the tool axis = ^A + A + A 2  vector  in Eq. (3.44) and its derivatives with respect to the spline  2  x  Z  parameter  must be computed. The resultant 1 , 2 st  nd  and 3 derivatives of rd  the orientation spline are expressed as:  dO  (3.46)  du  o. A  dQ 2  Oy  du  2  AA  o.  -2AA uu  u  AA ' "  Ul  A  3  °y  du  _  uuu\  3  2\A\A\ A\A\  I  \  uu\  u  I  ||  (3.47)  w_  Ul  UUU  3  _j_ 2  2  A \A\-A\A\  dQ  u  U  I  +2(A\A\  \m,  + A\A\  V »l  \uu  Iwl  \uu A  )A\ - 4A\A\ 2  »»l  \uft  I  \  "I  \A\)  2  \u I  1/  o. A\A\ "I  +A\A\ \uu  I  +2(A\A\ +2A\A\\A\  )A\ -6(A\A\ \A\ )  2  \uuu  V  " I  \u  3  I  3  \  I  I  IH  I  2  I  /  (3.48) As shown in Eq. (3.48), the third derivative of the angular spline in each direction depends on \A\, \A\ , \A\  and Ul  , which are functions of the upper and lower b-spline axis  components i.e. the control points ([P  X  P  P ] and [Q r  Y  z  X  Q  Y  Q^)-  Modifying any  control point of the upper b-spline in one direction affects the three axis components (O ,O ,O ) x  y  z  of the orientation spline since norm of tool axis vector (j^4|) changes. Therefore,  the problem is formulated as integral sum of all three derivative components of the orientation spline O ,O ,O xum  yuuu  zuuu  on the unit sphere. The new cost function for the angular  jerk minimization process can be written as:  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining Jo„ = p L u  +0  2  +0 p\u:  (3.49)  2  yuuu  49  zum  u  After the whole toolpath is fitted, the resultant knot vector is known and fixed. Control points of the upper spline ([Q  Q  x  Q ] ) are the free optimization parameters that can be T  y  z  used to modify the upper b-spline in order to penalize the cost function defined in Eq. (3.49) for the angular jerk. The vector of optimization variables is presented below;  0,  = [a a, - a. o  e,.  e„  '• Q  ya  Q  H  a, -  0.50)  Q J -  The constraints of the optimization problem are: (1) The max angular error, e  and  e  (2) the Euclidian distance between upper and lower b-splines {\A\).  The constraints on orientation error ef) are imposed in order to respect the machining (  tolerances and restrict the change in the upper b-spline. The limit on  will help to constrain  the relative distance between the splines near to H, which was pre-defined by the machining process. Since it is difficult to satisfy an equality constraint in nonlinear optimization process, a relaxation band is defined for \A\ of e = 0.01 ~ 0.25 mm so that the convergence is obtained faster by imposing inequality constraints instead of equalities. Given for e  e  and s, the constraints to be imposed are:  e - e  e  e  max  H-£<\A\<H  | + s\  Above constraints can be separated into individual inequalities;  (3.51)  , H  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining C (u ) e  = e -e  k  g  <0,  9max  C,-, (u ) = H-s-\A\  < 0,  k  (u ) = k  \A\-H-s<0  50  (3.52)  and they are grouped as a compact constraint vector as.  C{u)=[c  C,-, C{J  eo  for0<i/<l.  (3.53)  The optimization problem is then (Eq. (3.49) and Eq. (3.53))  J  n  =  min  ~\o  U  0  2  0 Yu  2  xmu+  Subjectto: C( )<0 for 0 < u < 1 (3.54)  2  yuuu+  zuuu  W  [Q, Qy &]u=0  Since the objective function and the constraints have nonlinear relationship with the design variables, the formulation in Eq. (3.54) results in a nonlinear optimization problem unlike in Eq. (3.36). It is solved using M A T L A B ' s Optimization Toolbox. The cost function and its gradients can be calculated analytically when the basis functions are given in analytical form after b-spline fitting for the angular optimization. Otherwise, numerical integration methods are used to evaluate the integral in Eq. (3.54) where adaptive Gauss Quadrature or a simpler method such as Simpson's 3/8 can be employed. Constraints are checked only at pre-calculated parameter values (w) obtained from Eq. (3.28). The gradient of the constraints can easily be supplied analytically to the algorithm, since they are simply a function of basis functions, N.Ju),  . For the algorithm, the evaluated basis functions are  •P \u=u  k  computed once and stored to be called by the algorithm each iteration.  3.7. Simulation Results  The feasibility of the toolpath generation algorithm presented in this thesis is shown using simulations. The aim of the first step is to fit the double b-spline toolpath to given discrete toolpath points (G01 Commands) with respect to specified tolerance limits. The  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  51  machining tolerance applied on aeronautical parts is given around ~0.025 [mm] [1]. The conventional spline fitting algorithms build in CNCs for high speed machining recommend fitting tolerances 0.1 [mm] for rough machining and from 0.05 [mm] to 0.01 [mm] for finishing operations of free form surfaces [69]. In order to show the potential of the presented fitting algorithm, two fitting tests are presented in the following. The b-splines of the generated toolpaths are selected at 5 degree. th  The error limit for the first b-spline toolpath generation test is set 0.02 [mm] for tool tip and 0.001 [rad] for the tool axis orientation. The test toolpath is initially discretized by C A M software for five-axis machining using 584 G01 points. The discrete points of the toolpath and the fitted b-splines can be seen on the Fig. 3.14 (a). Since the error bounds are selected tight, the curve and points are overlapping. The fitting process started with 29 initial control points and required -2.82 [min] on a P4 2.6 desktop computer running M A T L A B to result a b-spline toolpath consisting of 42 control points. The error distributions for the tool - tip positioning and tool axis orientation are given in Fig. 3.15 (a). On straight sections, the fitting error for the tool - tip drops down to 0.002 [mm] and less than 0.0001 [rad] for the tool orientation, which indicates that especially on smooth parts, tighter tolerances can be achieved. However, towards the end of the toolpath where there is sharp change in the curvature, the fitting errors are on limits. The resultant maximum positioning error is found max 0.018941 [mm]  as  e  =  for  the  tool  tip  and  orientation  error  is calculated  a s e ^ = 5.22 x 10^ [rad] = 0.029 [deg] using Eq. (3.32) and Eq. (3.35). The toolpath for the second test is shown in Fig. 3.14 (b). It is a round shape with constant change in the curvature that is initially defined with 183 TCP. The b-spline fitting algorithm shown in Fig. 3.15 started with 9 control points and required 25 control points to ensure specified 0.01 [mm] tool tip and 0.0005 [rad] tool orientation tolerances. Unlike the first fitting test, the errors are distributed evenly as shown in Fig. 3.15 (b). The maximum tool tip error is computed as e  max  e  = 4.58xlO" [rad] = 0.026 [deg]. 4  e  = 0.009 [mm] and the max tool orientation error is found as  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  52  The results indicate that very tight tolerances can be achieved by the proposed fitting algorithm in Section 3.5. Especially, for roughing and pre-finishing operations, the method allows generation of b-spline toolpaths for high speed five-axis machining.  Figure 3.14: Test Double B-spline Toolpaths Generated from Discrete Points  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining a) Fitting Errors of Tool-path #1  CO s  o s s  LU  CO l— _  t  03 L—  o LU  0.2  0.4 0.6 Spline Parameter  0.8  b) Fitting Errors of Tool-path #2  E CO  o  I— 1_  LU  CO  i—_ CO l—  o 1  LU  04  06  Spline Parameter Figure 3.15: Fitting Errors of the B-Spline Toolpaths  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  54  The effect of smoothing objective presented in Section 3.6.1 is investigated in the following simulations. Firstly, the influence of smoothing coefficients is studied on 45 two dimensional data points shown in Fig. 3.16 for an X - Y table. The fitted b-splines have 20 control points and constrained at both ends to the data points. In Fig. 3.16 (a) a regular least squares fit without any smoothening affect can be seen. The curve is deteriorated especially at the start section because of the incompetent knot distribution estimated by arc the length approximation. The average error of the fit is e  mean  =  0.211 [mm]. In Fig. 3.16 (b), the presented jerk minimization objective is  implemented into the least squares algorithm using Eq. (3.43). Fitting and smoothening coefficients are selected as y = le6 and a = 0.01 respectively. As a result, he curve shows higher average fitting error of e  mean  = 0.404 [mm]. But the oscillation at the beginning is  eliminated and it displays improved change in the curvature (See Table 3.1). Furthermore, it can be observed by eye-inspection that the optimized curve is smoother. In Fig. 3.16 (c), the fitting coefficient is increased to y = lei. As result of this, the average error is decreased to e  mean  =  0.255 [mm] with some loss of the curve smoothness. In the last simulation Fig. 3.16  (c), the fitting weight is further increased to y = leS, which reduces the average errors down to e  mean  =0.217 [mm]. The resultant b-spline is almost identical to the regular least squares  fit. Furthermore, in Fig. 3.17, the third derivatives of the fitted curves can be compared to judge the smoothness quality of the curve. Reducing the fitting coefficient diminishes clearly the fluctuations in the 3 derivative of the spline resulting in a smoother curve. As a result, rd  the acceleration components (2 derivative) on the b-spline ( See Fig. 3.17 (a) and Fig. 3.17 nd  (b) ) contains less change. As observed from the simulations, the most suitable coefficients can be selected as y = lei and a = 0.01 in order to smooth the curve while still preserving an acceptable fitting effect. The comparison of curvature changes is presented in Table 3.1.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining —  55  Optimized Least Squares Fit — Regular Least Square Fit » Shape Points Knot  c) gamma = 1e7, alpha = 0.01, Mean Error = 0.255[mm]  d) gamma = 1e8, alpha = 0.01, Mean Error = 0.217[mm]  Figure 3.16: Multi Objective Least Squares Smoothening of a B-spline - a) Regular Least Squares Fit, b-d) Optimized Least Squares Fit Comparisons  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  Integral of Rate of Change of the  Fitting & Smoothing  Curvature  Coefficients  ((/,«)) y = 1, a =  0  1.409 e  Y = le6,a = 0.01  8.57 e  y = \e7,a = 0.0\  1.259 e  y = \e8,a = 0.0\  1.04 e  17  15  16  17  Table 3.1: Comparison of Curvatures y=1e8,a=0.01  0  - - - y=1e7,a=0.01  — y=1e6,a=0.01  a) 3  0 >  > a3 •  c Q. W  C)  d)  3  3 "a  5?  x 10  Q)  CD  Q  Q  CD  CD  c  C  Q.  Q.  CO  02  0.4  0.6  0.8  Spline Parameter  02  0.4  0.6  0.8  Spline Parameter  Figure 3.17: 2 and 3 Derivative Comparison of Smoothened B-splines nd  rd  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining .  57  The smoothing effect is applied to the generation of lower b-spline for the first test toolpath shown in Fig. 3.14 (a). The weights for fitting and smoothening effects are chosen as y = \e7 and a = 0.01 while setting the maximum error tolerance as e  max  = 0.02 [mm].  However, using the above fitting and smoothing coefficients, the algorithm did not converge to satisfy the specified tolerance band. Consequently, the weight for smoothing effect is further reduced to or = 0.001. The smoothness of the b-spline fits for the toolpath is measured by comparing the 3 derivatives of regular and optimized fits shown in Fig. 3.18. rd  As inspected, the quality of fitted b-splines is very close. Because of the tight error limitations and reduced smoothening coefficient, both of the b-splines follow almost the same pattern. As a result, tuning of the smoothening effect is quite difficult in case of tight tolerance specifications and larger number of tool center points. For machining operations, the fitted b-spline is overly constrained by tolerance requirements allowing less change in the quality. Therefore, the optimization technique should be employed for generation of toolpaths with relaxed tolerance specifications.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  >  CO  0 Q 0  x 10  2  CD  Optimzed - -  5  58  Un-Optimzed  o T3 X %  -2  "5.  x 10  0 >  0.1  02  0.3  0.4  0.5  .0.6  0.7  0.8  0.9  0.1  02  0.3  0.4  0.5  0.6  0.7  0.8  0.9  0.1  02  03  0.4  0.5  0.6  0.7  0.8  OS  5  ro 0 "O Q >, 0^ "Q.  CO  0 >  -10  0  x 10 4r  5  ro >m 0 "O Q "N 0<^  "5. CO  Spline Parameter [u] Figure 3.18: Comparison of Optimized Lower B-Spline  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  59  The angular jerk optimization technique presented in. Section 3.6.2 was applied on the first b-spline toolpath shown in Fig. 3.14 (a) to inspect its effectiveness. Optimizing this long toolpath at once proved to be a very slow and inefficient process when considering the large number of control variables and constraints. In order to reduce the total processing time of the non-linear optimization, a technique similar to windowing scheme was implemented. The local modification property [58] of b-splines makes it possible that modifying a control point Pi affects the curve only in the parameter range of u e [p.  p  i+p+]  ]. The upper b-spline is  optimized locally at each step without altering shape of the whole curve. Hence, an optimization window size is set as with 5 control points. The window is shifted from the start towards the end of the toolpath while optimizing the section affected only by those control points. In case of non-convergence; the size of the window was increased by including one more control point from the un-optimized region. The optimization procedure for this toolpath consists of 42 variables described above and 3x584 constraints to be checked. Using the windowing technique; ~5.1 [mins] are required on P4 2.6 desktop computer running M A T L A B resulting in a feasible offline optimization procedure. The results for the optimized kinematic profiles of the tool at constant feedrate of 35 [mm/sec] are presented in Fig. 3.19. The angular acceleration profile presents fluent changes with less fluctuations resulting in a smoother trajectory. Hence, less vibration components are injected into the system while improving the overall tracking performance. The mean angular acceleration and jerk values at constant feed are smoothened and reduced by 9 [%] and 37 [%] respectively. The max angular acceleration is reduced by 54 [%] and the max angular jerk is brought down 79 [%] after the optimization. The effect of the angular jerk optimization on individual rotary axis profiles is presented in Fig. 3.20 and smoothened linear axis profiles are given in Fig. 3.21. Spikes in the linear axis acceleration and jerk profiles causing possible actuator saturations as well as excessive tracking errors are annihilated. For instance, the highest reduction is observed in the Y axis jerk profile by 80 [%]. Results indicate that the angular jerk optimization process reduces the overall jerk content of the toolpath efficiently allowing operating the machine tool at higher feedrate without violating the drive's kinematic limits.  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  40 I 0  i  i  i  i  i  1  2  3  4  5  Time[sec]  Figure 3.19: Angular Kinematic Profiles  I  6  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  a) Acceleration A Axis  'b  cu .52  1 0.5  r  Un-Optimized Optimized  0  E -0.51-1  b) Jerk A Axis o cu  20 10 0  CD  -10 -20  c) Acceleration C Axis 2 1 0 -1 -2 -3  d) Jerk C Axis o  40 20  lo -40 -60 2  0  s  Time[sec]  Figure 3.20: Rotary Axis Kinematic Profiles  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  Figure 3.21: Linear Axis Kinematic Profiles  62  Chapter 3. Spline Toolpath Generation from Discrete Toolpath Data for 5-Axis Machining  63  3.8. Summary  A continuous smooth representation of the toolpath is significant for higher part quality in high speed five-axis machining. In this chapter, a continuous toolpath consisting of double b-splines is generated from discrete tool path data. The fitting errors are limited depending on the given number of points and toolpath complexity as low in the range of 0.01 - 0.05 [mm]. Therefore, the implementation of such offline b-spline toolpath generation is suitable for most rough, as well as various finishing operations. In addition, smoothing of the toolpath is investigated for part quality and tracking performance by utilizing jerk objectives. The first optimization scheme for the tool tip spline did not found feasible because of tight tolerance specifications applied on the curve. However, the angular jerk of the toolpath is reduced using nonlinear optimization procedure. It is applied after the toolpath is fit and can handle long toolpaths time efficiently. The resultant angular acceleration profile shows smoother trends that enhances the tracking performance of the machine tool, allows increasing the overall feedrate.  64  Chapter 4 Interpolation and Feedrate Profiling for 5-Axis Machining  4.1. Introduction  A toolpath in terms of b-splines is defined for five-axis machine tools in Chapter 3 of the thesis. The next step is to generate an optimum feedrate profile for five-axis machine tools processing this b-spline toolpath, and interpolate those feed commands to the machine tool axes. The real time interpolation of b-spline toolpaths with minimal feed fluctuation is investigated in this Chapter. Furthermore, a jerk continuous and smooth cycle time optimal feedrate profile is generated for five-axis machining. Feedrate fluctuations arising from incorrect spline parameter calculation is discussed in the Section 4.2, and the 'Feed Correction Polynomial' method, developed by Erkorkmaz and Altintas [22] is implemented on b-spline toolpaths for feed consistent interpolation. Results from the 'Feedrate Correction Polynomial' approach are compared to the widely used, 'Taylor Series Expansion' [50] method in the simulations section. The Taylor approximation is also briefly explained. A feedrate scheduling algorithm was proposed by the same authors processing quintic spline toolpaths on Cartesian machines [3]. In the feedrate profiling work, the algorithm is modified for the double B-Spline toolpaths and implemented considering the 5-Axis machine tool's kinematic transformation and drive limits. The algorithm generates smooth kinematic profiles that consider machine tool's kinematic limits such as velocity, acceleration, jerk for minimum cycle time while aiming for better tracking performance.  4.2. Interpolation with Minimal Feed Fluctuation  For a specified feedrate profile V = V(t), a resulting displacement profile s = s(t) is obtained. Then, the position commands to the drives are calculated in order to track the  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  65  specified displacement profile during the machining operation. This task is performed by the real time feedrate interpolator embedded within the machine tool controller. The displacement profile can be provided in a linear, circular or a spline toolpath [42] format. For any of these provided toolpaths, the feedrate interpolator must calculate the axis increments as reference commands to the control loop at each sampling interval (j^). This task becomes computationally expensive when processing spline toolpaths [74]. In Fig. 4.1, a parametric b-spline toolpath is illustrated with a given displacement profile. As the spline parameter is incremented (Au) this corresponds to an arc increment (As) on the b-spline. The travel distances (i.e., the arc increments, (As)) at sampling intervals are scheduled during the feedrate generation step. Hence, the task of the real-time interpolator is basically to increment the spline parameter u = u(t) accurately to satisfy the desired arc increment s = s{t).  Figure 4.1: Real-time Spline Interpolation  In particular, b-spline toolpaths are defined in the parameter interval of we[0  l]. A  straightforward approach for real-time interpolators is to linearly increment the b-spline parameter proportional to the ratio of the arc increment to the total arc length (S ) (i.e., lotal  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  Au -l/S  66  • As). However, this linear relationship between spline parameter (w) and arc  l0!al  displacement (s) is not valid over the entire b-spline [21]. An imperfect calculation of the spline parameter would cause fluctuations on the actual feedrate. Such fluctuations would result in higher accelerations and jerk components on the drives, causing unwanted vibrations of the tool and deterioration from the part quality.  4.2.1. Taylor Series Approximation  The Taylor series approximation for real time spline interpolation was first introduced by Shipitalni et. al. [67] to estimate the next point along a spline toolpath on constant feedrate. In the literature, this method is widely used for interpolation of the b-spline toolpaths [74, 87, 92, 93]. It is based on calculating the curve parameter by using the differential relationship between the arc length (s) and spline parameter (w) when the sampling interval is sufficiently small. In the practice of five-axis machining, the feedrate is practically defined on the tool tip, along the lower b-spline (tool tip spline) of the toolpath P(u), the feedrate can be related to the spline parameter (w) using differential geometry as follows:  v(th  ds _ ds dw  du _  dt  dt  du dt  fit)  ds  (4-1)  du  The differential relationship between the arc length and spline parameter on P(u) = [x (w) y' (w) z' (w)] can be expressed as: 1  ds = dw  dP{u) du  1  T  du = \P (w \ du = ^(x'Xu)) + {y' {u)) +(z' {u))du U  u  u  (4.2)  The subscript variables are used to denote differentiation with respect to that variable as,  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  x:=dx'/du  =  ZNl:l(u)-P =Z  'dtf,»  i=0  X = d y / d u = x^S(")-^  xt  du  Z  du  "  du  '  =  z;=dz'/du=5:^S(«)-^ = Z ;=0  i=0  >  p  /=0  xi  67  (4.3)  y i  z>  and similarly,  x:„=dv/du =x<)( ).p,=x 2  du 2  M  i=0  1=0  V  ATI  y  ^ 2  ^=dV/du =t<)(ii).^=£  (4.4)  2  1=0  z  1=0  :„=dv/du =x^ »-p,=l: 2  2  where N (u) ip  du v — AJ2 .A d w , / '» 2  2  y A  AT  du 2  '  zi  is the p degree basis functions, N\ \u) is the k derivative computed as th  k  th  p  expressed in Chapter 3 by,  M i-"  (4.5)  N, Mi+ +\ ~ Mi+i i+P+  MM'Mi  and  d^N («)  rd^-'w,^( )/d ^ ) [  Lp  dw  P  w  =P  Mi  +P  W  d^  -Mi  l + 1  M  ^M/W*-  i+P+i  (4.6)  ~Mi-  Substituting Eq. (4.2) into Eq. (4.1) yields,  dw dt  V(t)  M»)J+k(4'  (4.7)  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  68  A computationally efficient solution for Eq. (4.7) is necessary in order to achieve a realtime interpolation of the parametric curve; an approximation based on Taylor's expansion is used. For instance, when the curve has a large radius of curvature and the sampling interval is sufficiently small enough (e.g., T = 1 msec) a first order approximation is usually adequate s  at constant feedrate;  , rp du dt  V-T.  = U; +  (4.8)  However, when the curve has small radii of curvature, the second order approximation is necessary, which is expressed as,  T  2  +-  du 2  (4.9)  2 dt  2  where (d u/dt ) is calculated by taking the derivative of Eq. (4.2): 2  2  du dt  2  du  dt'  .(4.10)  vit) (< (u)J  +  (y' (u)) u  2  (u)f f  +  A<{u))  (y[{u))  2  +  2  + (z[{uf  Substituting Eq.(4.10) into Eq.(4.9) yields the second order Taylor series approximation:  u  M  = u, +  V_J\  (V • T f x' { )x' («,)+ s  2  u  Ui  uu  y' { )yj{«,) u  Ui  +  z[{u )z' {«,) t  kk)) +ti(«,)J PM)  m  2  2  (4.11)  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  69  At each time step, the spline parameter is calculated by the first order approximation using Eq.(4.8), or for the second order approximation by Eq.(4.11). Following that, the calculated spline parameter is simply used in order to compute the position and orientation values on the b-spline toolpath for 5-Axis machining as mentioned in Chapter 3 of the thesis. However, the Taylor series expansion is not reliable when the feedrate is varying and results in an imperfect curve parameter calculation for toolpaths having sharp changes in the curvature. Furthermore, its real time implementation is quite computationally expensive for b-spline toolpaths [12] because of the calculation of basis functions and their derivatives. For instance, the basis functions and their k derivatives are calculated recursively by Eq.(4.5) th  and Eq.(4.6).  a)  b)  N (u) liP  U - U; U  U  i -1 " ^  U  + P  I+P -  i  N (u) liP  P-1  U  U  - i 1 U  + P  i  + P  -1"  -(P-1) i ,  u  U  +  i  N  lP  /  N . (u) i p  N (u) li0  2  N  i + 1 0  \  u  +  N  (u)  p  /  i + 1 p  N  i 1, -1( )  . (u) 2  N  i + p  .  \  N  j + 2 p  1 Q  (u)  +  (k-1)  (k-1)  N, .i(u)  - i 1 U  + P  i,p-l( )  i  N  U  +  1  , P »  . (u) 2  N  i + p 0  (u)  N . (u) i p  k  N  | + 1p  . (u) k  N  | + M p  _ (u) k  N  i + p p  . (u) k  Figure 4.2: Basis Function Computation a) Recursive Calculation of Basis Functions — b) Recursive Calculation of Basis Function Derivatives As illustrated in Fig. 4.2, the recursive computation of the basis function N (u) requires jp  (2 - l ) function calls to Eq.(4.5) [12]. In addition, regarding to the recursive computation P  scheme, the calculation the k derivative of a basis function N. (u) th  p  function calls to Eq.(4.6) and 2 (2"' k  k  requires (2* - l )  - l ) function calls to Eq.(4.5). Hence, for the k  th  derivative of a basis function, total number of {l - l ) function calls to Eq. (4.5) and Eq. p  (4.6) are required. Using Eq.(4.3), for a b-spline consisting of n + 1 control points, the  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  70  computation of a derivative at one point on the b-spline demands (« + 1^2 -1) function call p  to Eq. (4.5) and Eq. (4.6). The computational load is quite significant when the degree of the b-spline and number of axis are increased. For example, the computation of a point on a cubic b-spline (p = 3) with 5 control points (n = 6), 49 number of function calls are required where a quintic b-spline with the same number of control points demands 217 function calls.  4.2.2. Feed Correction Polynomial Method  Due to the deficiencies imposed by the Taylor's approximation, the feed correction polynomial method is introduced by Erkorkmaz and Altintas [22] for real time spline interpolation. In this subsection, the feed correction polynomial interpolation method is explained and applied specifically to the b-spline toolpaths. The feedrate is defined on the lower spline of the toolpath, which is a p degree b-spline th  parameterized by a knot vector (U) consisting of n + p + 2 knots (pi). Subtracting the p + l repetitive knots at the beginning and the end of the knot vector, there are in total  n-p  internal knots in U. This yields to N = n-p + l internal knot spans in the b-spline s  representation. The knot spans of the b-spline and their parameter ranges are given as  k=Q  .C/ = [0...0,  ^  p,  k  =  n  ~  p ,...,p _ ,p ,;  p  p+l  k=\  and illustrated in following Fig. 4.3.  n  l  n  P  ^  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  p+1  Knot Vector  71  P+1  Figure 4.3: B-spline Knot Spans and Segmentation  Since the parameter ranges are fixed by the consecutive knots, each internal knot span can be assessed as a piecewise spline segment. Therefore, the whole b-spline representation is divided into n - p +1 curve segments (See Fig. 4.3). As given in Eq.(4.12), the parameter range for the k segment of the b-spline is presented as: th  Mk+p — u — M i k+P+  where k = 0,..., n-p  .  (4.12)  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  The arc length, S  k  (M^lOO)  for k  th  72  segment is calculated by dividing each knot span into  subdivisions proportional to its .parameter range and computing the  corresponding axis positions as:  < =zZ i. (M +z-Au yp N  e  P  k+P  k  xi  i=0  Au =•  +p  k  (4.13)  I = 0,1,2,... M^ i=0 i=0  The segment lengths are computed by summing up the arc increments as:  n-p  (4.14)  'total ~ Z S k^k  k=0  The non-linear relationship between the b-spline parameter (w) and the arc increment (As) was already expressed in Eq.(4.2). According to the feed correction polynomial method, this non-linear relationship can be expressed by another polynomial. The objective is to use the feed correction polynomial in order to calculate the spline parameter as accurately as possible and with less computational load during real time interpolation of b-splines. Following that, the computed spline parameter is used to extract the orientation and position values on the double b-spline toolpath.  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  73  a) Feed Correction Polynomial XActual Displacement- Feed Correction Polynomial  Arc Displacement [s] b) Fitted Polynomials Over Knot Spans Hn+1 CD CD <D  .E  a. E  CO CD  CD CD  0-  0  0  [s]  b  i  "p+k Arc Displacement of Knot Span(k) -  Figure 4.4: Feed Correction Polynomial for B-spline Toolpath  The feed correction polynomial is illustrated in Fig. 4.4 (a). For a b-spline toolpath, it is specifically defined between the arc increments (As) and corresponding parameter changes (Au) of individual knot spans (segments) as shown in Fig. 4.4 (b). The feed correction polynomial segments in Fig. 4.4 (b) are connected to each other by continuity conditions in order to retain a single polynomial for the whole b-spline. Using seventh order polynomial spline [19], the feed correction polynomial for the k segment of the b-spline is written as:  u = A[s  1  +B{s  6  +...H{ (4.15)  where M  k+P  0 < s < S,  p  k+p+i  and  th  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  74  The feed polynomial defined in Eq.(4.15) basically supplies an accurate estimate of the b-spline parameter u to achieve the desired arc increment within each knot span during the interpolation. The parameter range (s) is bounded by the arc length of knot span (Sk). Considering Eq.(4.13) and Eq.(4.14) for the k  th  segment, the parameter increment and  associated arc travels are obtained offline, while computing the arc lengths as:  0  M  k+P  M  + Aw*  k+P  s  + 2Aw, ^>s =  u - M  k+P  2  _M _ S  Mk+ +\  k  P  (4.16)  where As ( ~ ijfakj k  k,e-l )  X  k « --y'k.t  +  e = z2&skn t=\  s  z  In order to avoid any ill conditioning, the parameter ranges in the feed correction polynomial are normalized by the arc length of each knot span (S ), k  u = a[a  +b[a  1  M  k+P  ^  u  6  +...h  J k  ^ Mk i  (4.17)  +P+  and 0< <7<1  where cr = s/S  k  and that results the upper bound of the arc displacement as 1. The  normalized polynomial coefficients are written as a{ = S\A  f k  , b{ = S\B{ ,...,/?/=  H(. The  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  75  first and second derivatives of the spline parameter with respect to the arc displacement (u - du/ds and u = d u/ds ) 2  s  are also expressed in normalized formulation:  2  s  u =f s  = ^[7a{cr 6b{^ ... ds bi.  {]  6  +  du  1  ds  Si  2  2  +  +  g  (4.18)  [42a{o- +30b{a +... 6  + 2f/]  5  The feed correction polynomial is simply obtained by applying a least squares fit over the normalized arc displacements and the spline parameter values (Eq.(4.16)) at each knot span as: ' 0  Mk+p «* = M  \  k+P+  +  2 A u  k  0  0  . ..  f  j • .. 1 H o-\ • .. 1  =  (4.19)  •. 1 1  1  1  . .. 1  H .  Referring to the least squares method, the objective function to be minimized is:  J = ±fe  =  (4.20)  ±(u-00 y(u-<M ). f  f  The first and second derivatives of the spline parameter with respect to the arc length can be written as:  u, -  du  du z  x' x' + y' y' +z' z' u  uu  u  uu  u  (4.21)  Therefore, the first and second order constraints at the beginning and the end of the k' knot span are given as:  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  mu  ,M,  U, =  76  _x x +y y +z z u  uu  u  uu  u  w  (4.22) final  u uu + yuy  X  final  -  X  m  "  + Z  n w Z  kJ+kf+kJ  The continuity boundary constraints have to be matched by the feed correction polynomial in order to connect individual polynomials fitted for knot spans over the b-spline toolpath. Accordingly, the ranges for b-spline parameter and arc travels are defined as zero boundary conditions on the polynomial. Those zero boundary conditions for k knot span are th  given as:  cr = 0 <=> u = p  '  k+p  a =1^ u=  (4-23)  p  k+p+i  Using Eq.(4.17) and Eq.(4.23), the zero, first and second boundary conditions are written in matrix form as:  'a{ 0  0  0  0  0 0 0  0  0  0  0  0 0  0  0  0  0  1  1  1  1  7  6  5  4  f  bi  1 0 ci 0 2 0 0 di  sW:;  1 1 1 1 ei 3 2 1 0 fl  s  42 30 20 12 6 2 0 0 si  s uT k  1  k U  (4-24)  r  Ci,. final  H Hence, the feed correction polynomial is obtained by solving the linear quadratic optimization problem similar to the one in the Chapter 3 by introducing the vector of Lagrange multipliers, A = [A  A\ ... A\] . The resultant linear equation system is: T  {  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  0 00 +L A T  = 0 u\  T  T  f  T  T  L  L6  r  ~0 U  L~  ~0 0  T  f  the  A{=a{/Sl.  solution,  the  B{=b{/Sl  (4.25)  0  Solution of Eq. (4.25) gives the coefficients 6 Following  77  of the feed correction polynomial.  coefficients  are  de-normalized  as;  H[=h{.  4.2.3. Iterative Solution of the Spline Parameter  The feed correction polynomial is defined by Eq.(4.17) and the coefficients of the polynomial {A(,B{,...,  are obtained from Eq.(4.25). However, when the feed  H{)  correction polynomial is directly used in the real time interpolation; it introduces imperfection into the b-spline parameter calculation because of the residual errors in the least squares fit. Therefore, the feed correction polynomial interpolation method can be combined with an iterative method in order to solve for the exact b-spline parameter during real-time interpolation as shown in Fig. 4.5 (c). a) Taylor's Approximation s(t). u  V(t)T  =u + i+i  i  u(t), B-Spline Toolpath  s  x (u>y (u>z ( ) 2  2  2  Ui  b) Feed Correction Interpolation s(t)|  u(t)=A s + B s + ... + H f  7  f  6  f  u(t) B-Spline Toolpath  c) Feed Correction with Raphson Iteration Newton lu s(t) u.init Raphson test B-Spline u(t)=A s +B s ...H Toolpath Iteration f  7  f  6  +  As des  Figure 4.5:  ^  A s  test  A  test  ' test y  ' test z  Arc | Estimate  Real Time Interpolation Methods  prev prev prev  78  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  From the feed correction polynomial, the initial guess u  inil  = u for the b-spline parameter  is obtained. Following that, the desired arc displacement As  for the control interval is  des  calculated using backward differencing of the displacement profile s(t). The arc increment As  lesl  corresponding to the tested spline parameter u  is estimated by computing the  test  resulting axis increments:  where: ^test. ~ test ~ prev • 4V'lest ' }'lest ~ }'prev X  X  a n t  ^ ^  test ~ test ~ prev Z  Z  n X'es, =Y. 'A >est) *i N  U  (4.26)  •  P  n  y'les,=Tu iA les,) yii=0 N  U  P  1=0  where x' , y' prev  prev  and z'  prev  are the axis commands applied in the previous control sample.  The error between the desired and tested arc increments is: (4.27)  and its gradient with respect to the tested spline parameter is derived using Eq. (4.26) and Eq. (4.27),  de  =  du,„<,  Ax' • (dxl, /du)+ Ayj, • [dy' /du)+ Az' , • Wtes, /du) lesl  tesl  As,„,  tes  (4.28)  79  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  where \dx' / du), \dy\ Jdu) and \dz' Jdu) are simply obtained by substituting the lest  calculated value of u  es  lesl  tes  into the derivative expressions given in Eq.(4.3) and (4.4). The  correct value of the spline parameter is computed regards the convergence of the below iteration (z):  k , , L =(»,«, ),•-/,  w  (4-29)  The desired convergence is achieved when the arc increment error e is below 10~ as in 8  Eq. (3.28). The fact that a very close guess («) is given by the feed correction polynomial and the analytical gradient (de/du ), the above iteration converges within 2 to 3 iterations, test  which is convenient for real time interpolation of b-spline toolpaths.  4.3. Feedrate Profiling for 5-Axis Machine Tools The feedrate profile of the trajectory is crucial for the quality of the machining process. In high speed machining, discontinuities and high frequency harmonics of the reference trajectory can easily excite the natural modes of the structure and servo system causing degradation from axis tracking and contouring performance. Therefore, the tangential feed profile defined together with the spline toolpath must ensure high continuity and smoothness allowing the axis servo drives to track the input reference kinematic profiles better. For basic linear and circular toolpaths, the feedrate profile consists of 3 main sections: acceleration, constant velocity (cruise) and deceleration. The starting acceleration and ending deceleration segments are defined with high smoothness, and their maximum values are determined with respect to individual servo drives acceleration and jerk limits. The cruise velocity is selected as high possible determined by the process parameters and axis kinematic limits aiming for the minimum machining time. On the other hand, b-spline toolpaths enable continuous change in the direction of the tool tip velocity and the tool axis orientation. Even interpolation at constant feedrate can produce excessive acceleration and jerk components on the drives when there is a sharp change in the tool orientation or in the spline curvature.  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  80  Therefore, the feed must be scheduled optimally in a smooth way throughout the spline toolpath to avoid violating any axis kinematic limits. For instance, the tangential feedrate at the tool tip must be reduced at sections where the change in the curvature is sudden to prevent exceeding acceleration limits of the drives and increased at straight sections to reduce the overall machining time. In addition, the non-linear inverse kinematics transformation from the workpiece coordinates to the axis commands must be considered for the five-axis machine tools. Minimum Jerk Double B-spline Toolpath  Velocity, Acceleration, and Jerk Limits  Toolpath Points  Feedrate Optimization  Solution of Quadratic Minimization Problem  Toolpath Description  Feed  Correction Polynomial 5-Axis Inverse Kinematics  I  P=P(u),0=0(u) Constant Arc Increment (As) Interpolation  x>  Minimum Cycle Time  itervt^... Minimum Jerk  Feed Profile  Axis Kinematic Profiles M=M(t),M,M,M Kinematic s, s, s  Profiles  Axis Geometric Profiles M=M(s), M , M , M s  sss  ssss  Figure 4.6: Overview of Feedrate Profiling In order to generate a continuously varying feedrate profile, Altintas and Erkorkmaz [3] defined the arc displacement profile as a jerk continuous (C ) fifth degree spline over the 3  specified travel distances on the spline toolpath, and optimized the feed profile for minimum cycle time. They also employed the minimum jerk criteria in the formulation for a smoother displacement profile. Similarly, a quintic displacement profile can be generated as a function of the time durations required to travel each knot span for the b-spline toolpaths. The overall feedrate profile can then be optimized by modulating the knot span durations in an iterative manner where the objective of the feed optimization is basically minimizing the total cycle time without violating any of the axis limits. The five-axis machine tool kinematic transformation is used to obtain the individual axis velocity, acceleration and jerk profiles on  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  81  the path. As a result of the iterative optimization process, a time optimal feedrate profile is generated that respects the velocity, acceleration and jerk limitations of the servo drives. The overall scheme of the feedrate profiling is shown in Fig. 4.6. 4.3.1.  M i n i m u m Jerk Feed Profile  As mentioned in the previous section, the complete b-spline can be partitioned into N =n-p s  + l parameter segments with arc lengths of S ,...S ...,S _,. 0  k  profile is then generated through the specified time durations,  The spline feed  N  ro,...,r ...,7^_, t )  required to  travel the corresponding knot span segment distances. The overall C displacement profile 3  fitted on the toolpath as a function time s = s(t), and the resulting kinematic profiles are shown in Fig 4.7. The displacement profile for the k knot span segment of the path is th  defined by a quintic polynomial of time as:  s {r) = A r  5  k  k  +B r  +C r  4  +Dr  3  k  +E r + F  2  k  k  k  ,0<r<T .  k  (4.30)  k  The resulting velocity ( i ) , tangential acceleration (s) and jerk ("i") profiles are obtained simply by differentiation with respect to r ;  s {r )=5A r  +4B r  4  k  S\(T)=  + 3C r + 2D r + E }  3  k  2  k  k  k  k  (4.31)  20 A r +\2B r +6C r + 2D 3  2  k  k  k  k  s (T) = 60A r + 245, r + 6C 2  k  k  k  In order to generate the jerk continuous (C ) feedrate profile, the piecewise fifth degree 3  displacement profiles for consecutive knot spans segments are connected to each other by imposing the position (s ), velocity ( i ) , acceleration (s) and jerk ('s ) boundary conditions,  s {T ) = s k  X0)\  k+  k  Sk fa ) =  5A T \ +4B T? +3C T  2  k  k x  k+]  +]  (°) ^ 20^ 7; + 12BJ 3  's'k fak) = ' ' i s  k+  (°)J  k  +2D T k  + 6C T + 2D  2  t  k  6 0 4 ^ + 24B T + 6C k  +i  k  k  k  k  k  +E =E k  k+l  = 2D  k+l  = 6C  k+l  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  82  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  83  Furthermore, as indicated in Chapter 3 of the thesis, smoothing of the generated profiles can be achieved by minimizing the overall jerk content. The objective function is the integral square of jerk throughout the whole motion expressed as:  (4.33)  Afc,-1  where T - £~^T is the total motion time. The overall jerk objective given in Eq.(4.33) z  k  is written for N local objectives corresponding to each segment as, s  Jv  = f D*o(j)] dr + ('[S;(r)] dr 2  2  + ...+ (*"[s^(r)f  dt  (4.34)  where the substitution of Eq.(4.31) in Eq.(4.34) results;  36004 r + 2%%0A B t + 2  Jv =] k  coefficients  (A , B ,C ,D ,E ,F ) k  k  k  k  k  k  3  k  (l20A C k  The  4  of  k  + 516Bl)u + 2S8B C r 2  k  the  k  minimum jerk  k  C  + 36C J  dr.  (4.35)  2  quintic  displacement  profile  segment are obtained by solving the linear quadratic optimization  problem similar to the Eq. (4.25), resulting from jerk minimization objective given in Eq.(4.35) and constrained by the continuity conditions (C°, C , C , C ) expressed by 1  2  3  Eq.(4.32). The complete formulation and generation of the minimum jerk displacement profile is presented in Appendix A of the thesis. The minimum jerk C displacement profile is obtained after the computation of the 3  coefficients. The objective of the feedrate optimization is then minimizing the total motion time (T ) within the kinematic limits of the machine tool. In order to optimize the total cycle E  time, a second cost function is formulated as the sum of the feed spline motion durations:  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  84  Ns-\  Tz = Z *-  (- )  T  4  36  4.3.2. Velocity, Acceleration and Jerk Constraints Constraints on the motion duration cost function in Eq.(4.36) are defined with respect to the kinematic limits of the five-axis machine tool's drives. For instance, the velocity constraints are defined to indicate the maximum allowable speed of each axis and the given tangential feedrate. In addition, in 5-Axis machine tool structures, the rotary axes have relatively low acceleration thresholds in contrast to Cartesian axes, which can cause actuator saturations at sharp curvatures and orientation changes. Therefore, acceleration limits are imposed to constrain the maximum reachable feed on the given toolpath. The jerk limits help to contain the frequency content of acceleration below the natural modes of the feed drives avoiding unwanted structural vibrations as well as excessive tracking error.  4.3.3. Inverse Kinematics of a Rotary Tilting Table 5-Axis Machine Tool The tool tip position [x  1  y*  z' f  and tool orientation [o  x  O  y  O J vectors are z  obtained along the double b-spline toolpath with respect to the generated position profile, s = s(t). Following that, the structure dependent inverse kinematics transformation must be employed to generate the Cartesian and rotary axis reference commands from position and orientation commands in real time. Commonly, five-axis machine tools have a kinematic structure based on a 3-Axis Cartesian machine with addition of two rotary axes as shown in Fig.4.8. The two rotary axes accomplish the desired orientation and the Cartesian axes perform the linear motion of the tool to the desired location. The five-axis machine tool configuration to be modelled is illustrated in Fig.4.8. It has a rotary-tilting table placed on the top of sliding the x-y table. The x-y table shows intended motion only along the X and Y directions, leaving the spindle to move the tool up and down. By utilizing a rotary-tilting table, two rotational degrees of freedom are added to the structure. The notation for this configuration is as follows: the rotary axis indicated with A  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  H5  (tilting axis) introduces a rotational degree of freedom around the X axis. On the top of it, the C (rotary axis) axis is located that introduces another rotational degree of freedom around the Z axis (tool axis). Since the cutting tool is symmetric about the tool frame, a five-axis machine tool can achieve any demanded orientation.  Figure 4.8: Kinematic Chain of the Rotary Tilting Table 5-Axis Machine Tool  In order to define the inverse kinematics transformation of this structure, the decoupling method [41] is used, and rotary joint angles are derived by the aid of homogeneous transformation matrices [51]. Firstly, frames (coordinates systems) are assigned for machine tool's elements such as the table, spindle, base and cutting tool shown (See Fig.4.8). The frame S is assigned for the sliding x-y table, which is connected to the base of the machine tool with a prismatic joint. The frame T is assigned for the tilting axis (A), which is a rotational joint fixed to the x-y table. Similarly, R is assigned for the rotary axis (C) attached to the tilting axis. W is assigned for the workpiece sitting on the rotary table. The frame SP  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  86  is assigned to the spindle that is connected to the column of the machine tool with a prismatic joint. It has a degree of freedom only along the Z direction. Another frame, OT is defined for the cutting tool. Finally, a fixed reference coordinate system, Ref is placed on the machine tool's base, which isrigidlyconnected to the column.  <—• X+5x  a)  b)  Figure 4.9: a) Translation Motion — b) Rotational Motion The rotational and translational relationship between these individual frames is presented using homogeneous transformation matrices. For example; pure translations along X, Y and Z directions at the amount of Sx, Sy, Sz between two frames is given by a homogeneous transformation matrix as,  H,translation  1 0 0 X + Sx 0 1 0 Y + dy 0 0 1 Z+& 0 0 0  (4.37)  1  where X, Y, Z are the initial offset values between the origins of frames. Considering  Eq.(4.37), the pure translation is illustrated in Fig.4.9 (a). The machining point is defined by the vector P = [x'  y'  z' ] with respect to the workpiece frame T  ) that is located on the  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  87  x-y table. Using Eq.(4.37), the transformed tool tip vector is represented with respect to the reference (Ref) coordinate system as:  F  y  =H  (4.38)  translation z'  1  Similarly, if the workpiece is located on the rotating table of the machine tool, the relationship between the reference and the workpiece frames is given by a pure rotation matrix. As shown in Fig 4.9 (b), a rotation around Z axis of the workpiece coordinate system by the angle of 86 can be represented by a pure 3x3 rotation matrix as, z  cos (6 +80 ) Z  -sin(d +86 )  Z  z  R = \- sin(6 + 50 ) z  z  0  .  0"  z  cos{0 +S6 ) z  0  Z  0  (4.39)  1  where 6 is the initial relative rotation between the frames. This rotation can also be defined Z  by the augmented 4x4 homogeneous transformation matrix as:  cos{e +5G ) -sin(9 +S9 )  0 0"  sin{d +S6 )  0 0  z  H.  z  z  z  z  z  cos[9 +S9 ) z  z  0  0  1 0  0  0  0 1  (4.40)  After the rotation is performed, the transformed tool tip vector in the reference coordinate system is expressed as:  88  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  x  P  y rotation  1  (4.41)  z  1  Respectively, intended rotations around X and Y axes are represented by similar homogeneous transformation matrices [51]. By the current frame approach, the homogeneous transformation matrices can be multiplied in a sequence in order to express the relative motion of the end frame with respect to any given reference frame. For instance, if N rigid bodies are connected in series and the corresponding homogenous transformation matrices between connecting axes are known, relation between the last coordinate system (N frame) th  in terms of the reference coordinate system is the sequential product of all of the homogenous transformation matrices:  N  H^Y\H^=H\HlHlH\...H  N-l N•  (4.42)  4.3.3.1 Rotational Motions The objective in five-axis machining is to keep the tool axis always aligned to the specified orientation and the tool tip placed at the reference machining point. In this modelled configuration, two rotary axes rotate the workpiece to align the orientation vector with the tool axis. Hence, the rotary axis rotation values are first calculated to establish the correct orientation of the tool.  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  + [0 0 1]  [O , O , O ] x  y  89  [0 0 1]  ^  z  . W  p  z  a)  Vef  x  b)  ref  Figure 4.10: a) Table Orientation Before Cutting — b) Table Orientation During Cutting In Fig. 4.10, the rotational motion of the table is shown. During the machining operation, the  workpiece  0(u)=[o (u) x  is  continuously  rotated  to  align  the  orientation  vector  O (u) O (u)f with respect to workpiece frame to the given tool orientation y  z  vector 0 '. The tool orientation is already fixed as 0 ' = [0 To0  To0  0 l ] along the Z axis with r  respect to the reference coordinate system. As shown in Fig. 4.10 (b), two rotations are performed in a sequence to accomplish this task. Thefirstrotation is carried out by the tilting axis A by the angle of 0 and followed by the rotary axis C by the amount of 0 . The A  C  rotation of the A axis is given by a homogeneous transformation matrix as:  1  0  0  cosii{0 )  0  sin(0j  0  0  A  0  0'  -sinfo)  0  cos(0j- 0 0  1  and the pure rotation of the C axis can be written as:  (4.43)  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining cos{0 )  -sin(0 )  0  0"  sin(0 )  cos(# )  0  0  0  0  1 0  0  0  0  c  c  c  c  H  e  =  c  90  (4.44)  1  Hence, the relative rotation between the workpiece and the reference frames is composed of these two sequential rotations. In order to achieve the correct orientation, the representation of the orientation vector respect to the reference coordinate system must equal to the fixed orientation of the tool vector. Using the current frame approach, the transformation matrices ( Eq. (4.43) and Eq. (4.44)) are multiplied in a sequence resulting in the following motion equation:  -H  H  0A  0  0C  =  0  T  O  0  (4.45)  ' .  Eq. (4.45) can be written as a set of equations in the matrix form as: -  1 0 0 0  0  0"  0  0 cos(e )-sin(6 ) sin{e ) cos(e ) 0 A  A  A  0  cos(e ) - sin[d ) sin(6 ) cos(e ) C  A  1  0  c  c  c  0  0"  0  0  0  0  1 0  0  0  0  1  "0" O {u)  QTOOI  0  1  1  y  0 {u) 2  1  (4.46)  1  Rotation angles 0 and 6 are obtained by solving the first two rows of Eq. (4.46), and A  C  this decoupling method gives solutions for the rotary axes shown in Table 4.1.  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  91  Table 4.1 Rotary Axis Joint Angles  o <o O <0 J0 (u)) tan ) ( +7t x  c  , tan  x  O >0  Undefined  Till  x  -It/2  _,(0l(u) 0l{u)) r^r +  y  O =0  x  e {u)  o >o  o,=o  y  tan~  [  {  O {u) z  x  tan  -A4  1  j0l(u) 0 (u)) tan -rx-— +  j  2  4.3.3.2 Cartesian Motions: Following the rotary joint angles, the Cartesian axis position commands are calculated to position the tool tip at the specified machining point. Homogeneous transformation matrices are defined to express the intended translations and the fixed offsets between the assigned frames. The time varying homogeneous transformation matrix between reference and x-y table frames can be written as,  "10 Ref  H  1 1  s  0  0 X(u) + S  x  1 0  0 0 1 0 0 0  Y(u)+S  y  (4.47)  S  z  1  where X(u) and 7(u) are time varying (w = u{f)) position of x-y table frame originfromthe reference frame origin. S , S , and S are the initial offset values between those frames, and x  y  z  can be selected depending on the initial position. Similarly, a fixed homogeneous transformation matrix is defined between table frame and tilting axis frame origins in terms of the offset values T , T , T : x  y  z  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  92  0 0 T  1  0 1 0 7 ;  (4.48)  0 0 1 7 ; 0  0 0  1  Another homogeneous transformation matrix between tilting and rotary axis frames is written as,  10 0 1 Hl = 0 0 0 0  where R , R , R x  y  0*, 0 ^ 1 R 0 1  (4.49)  z  are fixed offsets between frame origins. The workpiece is located on the  z  rotary table. The transformation matrix between rotary axis and workpiece frame origins is similarly given by:  1 0 0 0  Wp  H  0 1 0 0  0 Wp 0 Wp 1 W 0 1  x  (4.50)  y  Pz  where W , W and W are the offsets. The transformation matrix between the spindle and px  py  pz  tool frame is defined as:  TTSP  TJ OT nr 1  _  I  1 0 0 0  0 1 0 0  0 OT 0 OT 1 OT 0 1  y  (4-51)  z  where OT , OT and OT are the tool offsets. At last the homogeneous transformation matrix x  y  z  between spindle and the reference frame is given as:  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  Ref  H 1  1  SP  1  0 0  SP  0  1 0  SP,  0  0  1  Z(u)+SP  0  0  0  1  93  X  (4.52) z  where Z(u) is the varying spindle frame origin's distance from the reference coordinate system. SP , SP and SP are the initial offset values. X  y  Z  The resultant homogeneous transformation matrix between the fixed reference and varying spindle frames is written by the use of Eq.(4.51) and Eq.(4.52) as:  jrRef  OT  n  _  TrRef  TTSP  '  SP  (4.53)  OT '  Similarly, the homogeneous transformation matrix between the workpiece and reference frames can be written by employing the transformation matrices from Eq. (4.47) to Eq. (4.50), and Eq. (4.43), Eq. (4.44) in a sequence starting from reference to workpiece frames as:  w  H  f P  = s H  ef  -  H  T -  H  e  A  -  H  R -  H  t  (4.54)  JTR  The correct positioning of the tool tip to the given machining point is achieved by the motion of x-y table and the spindle. Consequently, the tool tip position is given in tool frame as P ° To  l  = [0  0  Of, and it is written with respect to the reference coordinate frame using  Eq.(4.53). The machining point P = [x'(w) y'(u)  z'(u)J  is also transformed from  workpiece coordinate system to the reference coordinate system using Eq.(4.54).  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  94  Frame Figure 4.11: Kinematic Transformation  Hence, the correct positioning is achieved as illustrated in Fig 4.11 by matching the representation of these two vectors with respect to the reference frame using Eq.(4.53) and Eq.(4.54) as:  "0"  x'{u) Ref  H  _  /(«) z'(u) 1  Machinig Position Vector  =  H  Ref  •  0 0  (4.55)  1 Tool Tip Position Vector  The solution of Eq (4.55) yields the resultant equations for the X, Y , Z Cartesian axis commands:  95  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining + sin{0M)-y'(")-^M-Wp  X{u)=OT +SP -R -T -S -cos{eM^'(u) X  X  +  Y{u) =  X  X  '  x  X  sin{0 {u))Wp c  0 T  y  + SP  Y  - T  Y  - S  Y  cos{e («)) • sin(6 (u)) • x' (u)-cos{0 («)) • cos{0 (M)) • y' («) +  -  Y  a  a  C  sin{0 (u))-x' (u)- cos(O («)) • sin{O («))• ^ a  a  + swfo  Z(u) = - 0 T  -  Z  • *Fp z  SP  + T  Z  + S  Z  cos(9  Z  +  A  ^  c  « / i ( 0 («)) • tf  z  fl  (u)) • sin(0 (u))-x'(u) + sin(6  sin{6  +  - COJ(^„ («)) • cos(O («))•  C  («))•  A  c  C  A  («)) • cos{d («))• / («) c  (u))-z' (u) + sin{e ( w ) ) • sin(6 (u)) • Wp + sin{O («)) • «w(0 («)) •  + cos(6  a  A  + cos{O («)) • ^ a  z  x  C  + s«(<9 («)) • i?^ + cos(6 0  A  c  A  ^  («)) • i?  z  (4.56) 4.3.4. Formulation of the Constraints The Joint variable vector,  M(U)=\X{U)  Y(ii)  Z(U)  6 {U) A  6 {U)J C  is obtained after  the inverse kinematics transformation using Eq.(4.46) and Eq.(4.56). Following that, the axis velocity, acceleration, and jerk profiles (M(t), M(t), M(t)) with respect to the displacement profile 5 = s(t) are expressed as,  M(t)  =  M(t)  =  M -s(t) s  M -s(tf+M -s(t) ss  M(t) = M -sit)'  +3M  sss  where M M , M s  ss  sss  (4-57)  s  SS  •s{t)-s(t)yM -sit) s  are the geometric derivatives of the toolpath computed over the axis  reference commands by  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  AM  M,(s) =  AM du  ds  dw ds d Mfdw^  dM  2  2  ds  2  dw  ds  dM  du  dw  \ds  2  dM 3  3  ds  3  3  For  96  2  du 2  + M.  (4.58)  ds  2  dw 3  the optimization process, the geometric derivatives of the toolpath are calculated  once at the initialization stage by interpolating the toolpath with constant arc increments (As), stored in a look up table and retrieved during the optimization process for faster evaluation (See Fig.4.6). The limits for tangential feedrate ^x ^y ,v ,v max  max  Zmax  ,v  0Amax  ),  acceleration  ffcmax  { j x m a x J Y r n a . J z m a x J e ^ a t J e ^ a * ) ,  a  r  e t  h  (a ,a ,a ,a Xmax  Ymax  Zmax  /  ,a  ()Amax  m a x  ),  , axis velocity and jerk  0cmax  nonlinear constraints to be imposed during the  e  optimization process expressed as:  0 < s(t) < f J max  ~  Xmax  V  ~  Ymax  V  — v, Z  max  ,  <Y(t)<v  ,  -  Xmax  £ Z(0 * v  ,  C  max  <e (t)<v  0<t<T . £  inequalities,  ecmax  ~  ^Z(t)<a ,  -a  ^' 8  c  rmax  eAmax  ^o (t)<v  a  £Y(t)<a ,  ~ Zmax  max  ~  Xmax  (t)  X  a  Zmax  A  A  "Xmax  - Ymax  Ymax  ^ 9 !  for  <X(t)<v  a  <0 (t)<a A  6A  JYmax  J Z max  Zmax  ffAmax  J Xmax  ~  — Y(t) ~  X(t) ^  ^  JYmax  ^(0 — JZ  max'  J 0 max  ~  max'  Jd  ~ ^c(0  A  c  max  Jxmax  ^A(0  (4.59)  max  —  —  j9 max A  J0  C  max  The above constraints are normalized and separated into individual  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  C}(0 =s(t)lf -\<Q,  C- (t) = -s(t)/f <0  max  c-(0 = X(t)/v  f  max  -\<0,C (t) = - X(t)/v  x  max  xmax  x  -1 < 0  Y(0/v c;(t)=  -1 < 0,  C-(t) = -Y(t)/v  Z(t)/v c (0 =  -1 < 0,  C~ (t) = -Z(t)/v  rim  +  zmax  z  • e (t)/ c;jo=  C y(0 = Y(t)/a +  C (t) = Z(t)/a z  C:  x  Xmax  C y(0 =Y(0/j -l< +  Z(t)la  c:  =-<? r  c  _  c  -1<o  -1^0 - i <o  =-x(t)ij  Xma  C ~ . f O = -Z(t)/a  0  -1 s o.  = w / ^  _ -1 * o,  c^/r;=- ^ r q  c  -  r o = - 0 ro/v, C  Jerk Constraints  -1 < 0  Zmax  „  Acceleration,  eAmax  rma  =^ ^/S,  -1 < 0  -1 < 0  C~f<; = -Y(t)/j -l<  -1 < 0,  Zmm  zmflI  A  - l <o, 0,  rma  = - Z(f ;/a  Sk  - 1 * o,  c (t) = x(t)/j  C O  Ymax  - I< o  -1 < 0  Xmax  c- (t) = -e (t)/v  max  cmca  0cmax  c  y  = $A(t)K -1<o, =0c(t)he  c  C~ (t) = -Y(t)/a  -1 < 0,  zmax  c- (t) = -e (t)/v x  -1 < 0,  Ymax  eAmax  A  C~ (t) = -X(t)/a  -1 < 0,  Xmax  dA  Velocity, -1 < 0  -I<o  c- (t) = -e (t)/v  - I < o,  V0c max  C\(t) = X(t)la  zmax  i  eAmax  c  -1 < 0  Ymax  -I<o,  e (t)/v  :  A  97  c mai  -  1  (4.60) *  0  -1 < o  The inequality constrains in Eq. (4.60) are grouped together and written in the vector form of  c(/)=[c;(0  c (t)... x  c (t)... x  ct(t)... i  c;(/)  c^)...f,(5.6i)  and formulated for the optimization process as:  C(/)<0  4.3.5.  for  o</<r }. 2  Solution of the Optimization Problem  The optimization problem is given by Eq. (4.36) and Eq. (4.62):  (4.62)  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  98  Ns-l  min T = min ]T T subject to C(t) < 0 for 0 < t < T , £  k  z  (4.63)  which is a linear programming problem with non-linear constraints because of the fact that the derivative of the cost function respect to the optimization variables are linear. It is solved by MATLAB's Optimization toolbox. Each segment duration is divided into (-5...100) subdivisions proportional to its corresponding arc length S  k  as shown in Fig. 4.12 and the  look-up table for the geometric derivatives is prepared by Eq. (4.58).  f Constraint Check  Figure 4.12: Constraint Evaluation  During the optimization process, the axis kinematic profiles are calculated using Eq. (4.57) over the look-up table and the constraints are checked at fixed subdivisions using Eq. (4.0) of the segment durations. As the segment durations are modulated, the number of constraints to be checked remains the same, but the evaluation locations change since the subdivisions are stretched out or compressed. Such a 'floating constraint' scheme allows analytical computation of the kinematic profile gradients with respect to the optimization variables (ds/dX, (dC(t)/dX),  ds/dX,  ds/dX,  ds/dX),  as well as the constraint gradient  which enhances the convergence of the optimization given in Appendix B.  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  99  4.4. Simulation Results  4.4.1. Comparison of the Interpolation Methods Four real time interpolation techniques of the b-spline toolpath are compared. The feedrate consistency and fluctuations are observed for the following interpolation methods:  (1) 1 order Taylor's Approximation st  (2) 2  nd  order Taylor's Approximation  (3) Feed correction polynomial interpolation (4) Feed correction with iterative interpolation  In the literature, the interpolation of b-spline toolpaths is performed using Taylor series approximations [74, 87, 92, 93]. Therefore, the performance of the feed correction polynomial interpolation technique is compared to 1 and 2 st  nd  order Taylor's Approximation  methods and summarized in Table 4.2. The b-spline toolpath for the tests is shown in Fig. 4.13. It is a round arbitrary designed shape with sharp changes in the curvature near quadrants. The toolpath is run at constant feedrate of 100 [mm/sec] and interpolated at 1 [msec] sampling. Results are summarized for the velocity of the tool tip in Fig. 4.14. Interpolation result of the toolpath with 1 order Taylor approximation is shown in Fig. 4.14 (a). The st  fluctuation on the tangential feedrate of the tool tip is observed as maximum 0.47 [%]. Since calculation of b-spline derivatives demand high computational load, each iteration needs ~ 4 [msec] on a P4 2.6 GHz desktop computer running MATLAB. However, it can be implemented by a much more efficient computational platform for real time interpolation. In Fig. 4.15 (b) results for the 2  nd  order Taylor Approximation can be seen. The maximum  fluctuation is reduced by 0.017 [%]. Since more terms are included in the approximation, the resultant time expense is increased to ~6 [msec]. In Fig. 4.14 (c), the same toolpath is interpolated using the feed correction polynomial. The correction polynomial is initially fitted offline and the real - time interpolation is performed without any recursive calculation of the spline parameter. Therefore, the computational load can be neglected (<0.02 [msec]). The maximum fluctuation in the feedrate is 0.03 [%], which is closer to the 2  nd  Order Taylor  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  100  approximation. This indicates that the method can be used for real time b-spline interpolation as compared to computationally expensive Taylor approximation techniques. In the last simulation, the toolpath is interpolated with the recursive interpolation combined with feed correction polynomial guess. There is practically no feedrate fluctuation for the b-spline toolpath. However, at least 2 iterations are needed, which increases the computational time up to ~8 [msec]. In contrast to any Taylor approximation method, it guarantees higher accuracy at for parameter calculation.  Table 4.2 Summary of Interpolation Results  [%]  Time Expense [msec]  1 order Taylor Series Apprx.  0.47  =4  2  0.017  =6  Feed Correction Poly.  0.03  <0.02  Feed Correction Poly, with iterative solution  0.0001  =8  Interpolation Method st  nd  order Taylor Series Apprx.  Feed Fluctuation  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  Upper B-spline  Figure 4.13: B-spline Toolpath  0)  LL  p £ ,  b) 2  Order Taylor Interpolation  c) Feed Coorection Polynomial Interpolation 100.04  To 8 100  »— Ui  •D  <u E «  £  02  100 99.98 99.96  d) Feed Coorection Polynomial Recursive Interpolation  0.3  0.4  0.3  0.6  0.7 .  e) Number of Iterations  0.8  0.9  1.0  1.1  1.2  1.3  Time[sec]  c o co cu 200  400  600  800  '  .1000  1200'  SampleQ  Figure 4.14: Feedrate Fluctuations  -.1400  1600  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  102  4.4.2. Comparison of Constant and Optimized Feedrate Profiles The feedrate optimization algorithm is tested on the toolpath that is generated in chapter 3 and shown in Fig. 4.15. It consists of b-splines defined by 42 control points which contains mildly varying curvature and orientation at the start section as well as sharp changes towards the end. The simulations are performed with respect to the kinematic limits of the HELLER MC 25 five-axis machining center. The velocity, acceleration and jerk limits are given in the following Table 4.3. The maximum feedrate limit on the toolpath is set as f  max  = 200 [ mm/sec]. The  optimization of the toolpath required 120 [sec] using a windowing technique [19] at P2.6 Ghz desktop computer running MATLAB and resulted in a motion duration of 2.09 [sec]. Trying to travel the same toolpath without violating any of the axis kinematic limits results in using a feedrate profile of 53 [mm/sec] with cubic acceleration/decelerations of 1250 [mm/sec ], and a tangential jerk of 15000 [mm/sec ] having an overall motion duration of 2  3  3.751 [sec]. With the optimized feed profile, the improvement on the total cycle time is 44 [%]. The feedrate profile is shown in Fig. 4.16, which indicates that the tool tip velocity can be increased up to 200 [mm/sec] at the start section of the toolpath while taking the advantage of the straight feasible region. On the other hand, it is scheduled conservatively towards the end where the change of tool orientation and curvature is quite sharp. The linear and rotary axis kinematic profiles are shown in Fig. 4.17 and Fig. 4.18 indicating that the major applied limits are the axis acceleration and jerk.  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  103  Table 4.3: Kinematic Limits of 5-Axis Machine Tool Heller MC25 Axis Limits Linear Axis(X,Y,Z)  Acceleration Limits  Velocity Limits v ,=1000  mm =10000—  W W  m  maz  Rotary Axis (A)  v  Rotary Axis (C)  v  e  e * cma  sec  sec . _ rad =sec 4  Amax  1  = 10  rad  2  n <:  r  a  5  sec  Jerk Limits  d  sec , _ . rad  «* «=17.3 cM  2  sec  u-  30000  sec  4ns - 62.8 sec - 9S1 Je max r  ^1  ^  ra<  3  sec  Figure 4.15: B-spline Toolpath  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  104  Optimized and Unoptimized Feedrate 200  Optimized Feedrate Profile Unoptimized Feedrate Profile /  2.09 Time[sec]  3.751  Figure 4.16: Optimized and Unoptimized Feedrate Profiles X AXIS  1000  Y AXIS  Z AXIS Limit Optimized,, .. . , \ Unoptimized  3.751  3.751 Limit Unoptimized  s ^•| V  1  2.09  Time[sec]  1  2.09  Time[sec]  3.751  I  Optimized 1 2.09 Time[sec]  Figure 4.17: Optimized and Unoptimized Cartesian Axis Profiles  Limit 3.751  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  A AXIS  C AXIS L/'m/r  4.2  «-—» o  105  jtimized  10.5 o  o  pptimized  TJmTf]  Unoptimized  Oj to  -4.2  Limit  1 |  6.3  V  CD O  o  3.751  2.09  Limit]  Optimized  -10.5  Limit]  0 17.3  Unoptimized  Limit. Optimized \  8n  3.751  2.09 unoptimized A  A  \  "5  CD O CO  o  L/'m/'f  -6.3 0 62.8 h  3.751  2.09 n  Optimized  ~LJmil \  Mm  Unoptimized ^ ^  -62.8  im  1  2.09 Time[sec]  3.751  -251  Unoptimized  Limit  1 2.09 Time[sec]  3.751  Figure 4.18: Optimized and Unoptimized Rotary Axis Profiles  4.4.3. Summary The interpolation of b-spline toolpaths and optimal feedrate profiling are investigated in this section. By utilizing the feedrate correction polynomial, it is shown that b-spline toolpaths are interpolated in real time efficiently. Feed fluctuations arising from inaccurate calculation of the spline parameter is further eliminated using the recursive interpolation scheme. For 5 - axis machining, the feed correction polynomial interpolation method allows tangential velocity at the tool tip to be kept consistent. However, the effective feedrate is changing along the tool axis because of the angular velocity of tool. Therefore, during flank milling operation, feedrate variations can be encountered along tool axis. In contrast, in 3 -  Chapter 4. Interpolation and Feedrate Profiling for 5-Axis Machining  106  axis machining this method can easily be employed for interpolation of any b-spline and Nurbs toolpaths in order to eliminate unwanted feed fluctuations. Furthermore, a feedrate scheduling algorithm for 5-Axis machine tool is implemented. The feedrate curve is defined as a fifth order C continuous spline providing jerk continuous 3  kinematic profiles to the drives. Furthermore, the optimized feedrate profile intends minimum cycle time while not violating the axes' velocity, acceleration and jerk limits. The algorithm can be used for high speed machining where continuously alteration of the feedrate profile is not crucial for the manufacturing process. For instance, its implementation is suitable for rough machining and pre-finishing of dies and molds combined with the previously presented b-spline toolpath.  707  Chapter 5 Conclusions 5.1. Conclusions In this thesis, a smooth and optimal trajectory generation algorithm is presented for fiveaxis high speed machining. A b-spline toolpath has been generated in order to deliver continuous tool motion for five-axis machining by compressing linear GI blocks assigned in the initial NC program. The real time interpolation of the b-spline toolpath is realized usigng the Feed Correction Polynomial method for minimal feedrate fluctuation. A jerk continuous feedrate optimization technique has been formulated for traveling along the b-spline toolpath in minimum time while adhering to the machine tool's kinematic limits such as velocity, acceleration and jerk. The thesis has the following contributions: The spline toolpath generation algorithm employs discrete toolpath position and orientation vectors extracted from the C A M system and generates a smooth toolpath composed of two b-splines with respect to the user given position and angular tolerance specifications. The double spline toolpath format is capable of interpolating tool position as well as the orientation vectors continuously in a synchronized way, decreasing the need for a machine tool dependent inverse kinematics post processing. The toolpath format is convenient for five-Axis high speed machining purposes and accepted by conventional CNC systems for real time interpolation. The smoothness of the generated toolpath is further improved using two optimization procedures. The first linear optimization technique is applied during the least squares fitting of the tool tip position spline, and it modifies the spline for smoother motion of the machine tool's Cartesian axes. However, the linear optimization is not practical when the tolerance specifications are tight. On the other hand, because of the machine tool's kinematic structure, any irregularity in the angular motion of the tool is directly reflected to the complete multiaxis motion. As an alternative the non-linear optimization algorithm addresses the angular motion of the tool, and eliminates any oscillation in the angular acceleration by modifying the orientation spline of the toolpath. The non-linear method the path to be more accurately  Chapter 5. Conclusions  108  tracked given the limited controller bandwidth, since any excitation of the servo and structural dynamics are avoided. The real time interpolation of the b-spline toolpath is improved using Feed Correction Polynomial Method. As compared to the conventional b-spline interpolation algorithms based on first order Taylor Series' expansion series; the implemented method is computationally less expensive and provides better feedrate consistency. Furthermore, combined with the iterative parameter correction scheme, the feed correction polynomial method exhibits an order of magnitude less feed fluctuation at steep curvature regions in contrast to second order Taylor's approximation. The feedrate optimization algorithm utilizes a quintic C spline as displacement profile in 3  order to realize a jerk continuous trajectory in five-axis machining. The spline displacement profile is especially constructed using minimum jerk criteria ensuring a smooth and continuous change in the feedrate enabling better tracking performance of the drives. The minimum cycle time objective is injected into the generation of optimal feed profile, and is constrained with the distinctive velocity, acceleration and jerk limits of the cartesian and rotational axes. The minimum cycle time method combines the ability to continuously modulate the feed, acceleration and jerk profiles along the spline toolpath for five-Axis high speed machining, and smoothness nature enhances the tracking performance of the servo control system while taking the advantage of machine tool's physical limits.  5.2. Future Research Directions The trajectory generation techniques presented in this thesis are practically applicable to the state of art five-Axis high speed machining. However, the proposed algorithms are tested only in virtual environment using the kinematic structures of widely used five-Axis machine tools. Experimental validations should to be conducted on physical test platforms, and the resulting contouring performance should be investigated including the tracking errors of each servo system using the nonlinear kinematic transformation between the joint space and workpiece coordinates. To expand the toolpath fitting technique to surface finish machining operations, limits for the position and orientation errors df the generated to.olpath must be readjusted by including  Chapter 5. Conclusions  109  the tool geometry, such as the tool radius of hemispherical or toroidal cutters in order to eliminate any possible interference between the tool and the workpiece. For instance, in the case of toroidal tool, the tool tip position error is coupled to the orientation error. The feedrate varies along the tool axis because of the angular velocity leading to fluctuations and oscillations in the feed vector during flank milling operations especially where the angular velocity is more dominant then the tool tip velocity. To improve the overall feed consistency, the upper b-spline of the toolpath must be re-parameterized or the interpolation step size is required to be adjusted that will reduce the angular acceleration of the tool. The feedrate profiling algorithm uses a smooth quintic polynomial spline. However, the replacement of quintic polynomial with a b-spline feed profile would enhance the convergence of the algorithm by utilizing its local modification property as an effective windowing technique. The actuation dynamics should also be included in feedrate optimization, which will prevent exceeding the actuation torques.  110  Bibliography [1]  Affouard, A., Due, E., Lartigue, C , Langeron, J.-M., Bourdet, P., 2004, "Avoiding 5Axis Singularities Using Tool Path Deformation", International Journal of Machine Tools and Manufacture, Vol. 44, No.4, pp. 415-425.  [2]  Alhanaty, M., Bercovier, M., 2001, "Curve and Surface Fitting and Design by Optimal Control Methods", Computer Aided Design, Vol. 33, No. 2, pp. 167-182.  [3]  Altintas,. Y., Erkorkmaz, K., 2003. "Feedrate Optimization for Spline Interpolation in High Speed Machine Tools", Annals of CIRP, Vol. 49, No. 1, pp. 265-270.  [4]  Altintas, Y., 2000, Manufacturing Automation: Metal Cutting Mechanics, Machine Tool Vibrations and CNC Design, Cambridge University Press.  [5]  Bedi, S., Quan, N., 1992, "Spline Interpolation Technique for NC Machines", Computers in Industry, Vol. 18, No. 3, pp. 307-313.  [6]  Bobrow, J. E., Dubowsky, S., Gibson, J. S., 1985, "Time-optimal Control of Robotic Manipulators Along Specified Paths", International Journal of Robotics Research, Vol 4, No. 3, pp. 3-17.  [7]  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.  [8]  Butler, J., Haack, B., Tomizuka, M . , 1998, "Reference Generation for High Speed Corrdinated Motion of a Two Axis System", Symposium of Robotics, ASME Winter Annual Meeting, Chicago, IL, USA, DSC Vol. 11, pp. 457-470.  [9]  Cao, B., Dodds, G.I., Irwin, G.W., 1997, " Constrained Time-Efficient and Smooth Cubic Spline Trajectory Generation for Industrial Robots", Proceedings of IEE, Control Theory and Applications, Vol. 144, No. 5, pp. 467-475  [10]  Cebula, A. J., Zsombor-Murray, Paul J., 2003, "Nurbs to Avoid Boundary Orientation Poses in Serial Manipulators", Journal of Robotic Systems, Vol. 20, No. 12, pp.723736  Bibliography [11]  111  Chang, Y.-H., Lee, T.-T., Liu, C.-H., 1988, " On-Line Cartesian Path Trajectory Planning For Robot Manipulators", Proceedings of 1988 IEEE  International  Conference on Robotics and Automation, pp. 62-67. [12]  Cheng, M.-Y., Tsai, M . - C , Kuo, J.-C, 2002, "Real-time Nurbs Command Generators for CNC Servo Controllers", International Journal of Machine Tools and Manufacture, Vol. 42, No. 7, pp. 801-813.  [13]  Chou, J.-J., Yang, D.C.H., 1991, "On the Motion Coordination and Planning of Multi Axis Manipulators with Orthogonal Regional  Structures",  Fifth  International  Conference on Advanced Robotics, pp. 1814-1819. [14]  Chou, J.-J., Yang, D.C.H., 1992, "Tool path Geometry for Multi-axis Kinematics Conversion in CAD-CAM Integration", ASME Design Engineering Division, 22nd Biennial Mechanisms Conference, Vol. 46, pp. 519-529.  [15]  Chou, J.-J., Yang, D.C.H., 1993, "Automatic Generation of Piecewise Constant Speed Motion for Multi-axis Machines", ASME Design Engineering Division, 22nd Biennial Mechanisms Conference, Vol. 46, pp. 531-536.  [16]  Constantinescu, D., Croft, E.A., 2000, " Smooth and Time-optimal Trajectory Planning for Industrial Manipulators Along Specified Paths", Journal of Robotic Systems, Vol. 17, No. 5, pp. 233-249  [17]  Elkeran, A., El-Baz, M.A., "Nurbs Feedrate Adaptation for 3-Axis CNC Machining", Maintence Resources, http://www.maintenanceresources.com / referencelibrary/ ezine /nurbs.htm  [18]  Erkorkmaz, K., 1999, M.A.Sc. Thesis:" High Speed Contouring Control for Machine Tool Drives", University of British Columbia, Department of Mechanical Engineering, BC, Canada  [19]  Erkorkmaz, K., 2004, Ph.D. Thesis:" Optimal Trajectory Generation and Precision Tracking Control for Multi-Axis Machines", University of British Columbia, Department of Mechanical Engineering, BC, Canada  [20]  Erkorkmaz, K., Altintas, Y., 2001, "High Speed CNC System Design: Part I - Jerk Limited Trajectory Generation and Quintic Spline Interpolation, " International Journal of Machine Tools and Manufacture, Vol. 41, No. 9, pp. 1323-1345  Bibliography [21]  112  Erkorkmaz, K., Altintas, Y., 2001, "Trajectory Generation for High Speed Milling of Molds and Dies", Proceedings of the 2nd International Conference and Exhibition on Design and Production of Does and Molds, Kusadasi, Turkey, D M 4 6 .  [22]  Erkorkmaz, K. Altintas, Y., 2002, "Optimal Spline Interpolation for High Speed Machine Tools", Proceedings of 6th Biennial Conference on Engineering Systems Design and Analysis, Istanbul, Turkey.  [23]  Erkorkmaz, K., Altintas, Y., 2003, "Quintic Spline interpolation with Minimal Feed Fluctuation", Proceedings of the ASME Manufacturing Engineering Division, MEDIMECE2003-42428.  [24]  Dyllong, E., Visioli, A., 2003, "Planning and Real time Modification of a Trajectory using Spline Techniques", Robotica, Vol. 21, pp.475-482.  [25]  Gossard, F.-L., David, C , 1995, "Multidimensional Curve Fitting to Unorganized Data Points by Nonlinear Minimization", Computer Aided Design, Vol. 27, No.l, pp. 48-58.  [26]  Farin, G., 1993, Curves and Surfaces for Computer Aided Geometric Design - A Practical Guide, Academic Press.  [27]  Farouki, R. T., Shah, S., 1996, "Real-time CNC Interpolators for Pythagoreanhodograph Curves", Computer Aided Geometric Design, Vol. 13, No. 7, pp. 583-600.  [28]  Flash, T., Hogan, N., 1985, "The Coordination of Arm Movements: an Experimentally Confirmed Mathematical Model", Journal of Neuroscience, Vol. 5, pp. 1688-1703.  [29]  Fleisig, R.V., Spence, A.D., 2001, "Constant Feed and Reduced Angular Acceleration Interpolation Algorithm for Multi-axis Machining", Computer Aided Design, Vol. 33, No.l, pp. 1-15.  [30]  Goldenthal, R., Bercovier, M . , 2004, "Spline Curve Approximation and Design by Optimal Control Over the Knots", Computing, Vol. 72, No. 1-2, pp.53-64.  [31]  Hohenberger, W., Reuding, T., 1995, "Smoothing Rational B-spline Curves Using the Weights in an Optimization Procedure", Computer Aided Geometric Design, Vol.12, No. 8, pp. 837-848  [32]  Horn, B. K. P., 1983, "Curve of Least Energy", A C M Transactions on Mathematical Software, Vol. 9, No.3, pp.441-460.  Bibliography [33]  113  Jee, S., Koo, T., 2003, "Toolpath Generation for NURBS Surface Machining", Proceedings of the American Control Conference, Vol. 3, pp. 2614-2619.  [34]  Jensen, C.G., Anderson, D . C , 1992, "Accurate Tool Placement and Orientation for Finish Surface Machining", ASME Production Engineering Division, Winter Annual Meeting in Anaheim, CA, USA, Vol. 59, pp. 127-145.  [35]  Jouaneh, M., Wang, Z.-X., Dornfeld, D., 1988, " Tracking of Sharp Corners Using a Robot and a Table Manipulator", Proceedings of the USA - Japan Symposium on Flexible Automation, pp. 271-278.  [36]  Juttler, B., 1994, "Visualization of Moving Objects Using Dual Quaternion Curves", Computers & Graphics , Vol. 18, No.3, pp. 315-326.  [37]  Kim, J.-H., Ryuh, B.-S., Pennock, G.R., 2001, "Development of a Trajectory Generation Method for a Five-axis NC Machine", Journal of Mechanism and Machine Theory, Vol. 36, No.8, pp. 983-996.  [38]  Koninckx, B., Van Brussel, H., 2002, "Real-time NURBS interpolator for distributed motion control", Annals of CIRP, Vol. 51, No.l, pp. 315-318.  [39]  Koch, P. E., Wang, K., 1988, "Introduction of B-splines to Trajectory Planning for Robot Manipulators", Modeling, Identification and Control, Vol. 9, No.2, pp. 69-80.  [40]  Koren, Y. Lin, R. S.,1994, "Real - Time Five Axis Interpolator for Machining Ruled Surfaces", Proceedings of the 1994 International Mechanical Engineering Congress and Exposition, ASME Dynamic System and Control Division, Vol. 55-2, pp. 951959.  [41]  Koren, Y.,1995, "Five Axis Surface Interpolators", Annals of CIRP, Vol. 44, No. 1, pp. 379-382.  [42]  Koren, Y., Lo, C C , Shpitalni, M . , 1993, "CNC interpolators: Algorithms and Analysis", ASME Production Engineering Division, Proceedings of the 1993 ASME Winter Annual Meeting, Vol. 64, pp. 83-92.  [43]  Kruth, J.-P., Kerstens, A., 1998, "Reverse Engineering Modeling of Free-Form Surfaces From Point Clouds Subject to Boundary Conditions", Journal of Materials Processing Technology, Vol. 76, No. 1-3, pp. 120-127  Bibliography [44]  114  Kyriakopoulos, K., Saridis, J., George, N., 1988, "Minimum Jerk Path Generation", Proceedings of the IEEE International Conference on Robotics and Automation, Philedelphia, PA, USA, pp.364-369.  [45]  Kyriakopoulos, K., Saridis, J., George, N., 1994, "Minimum Jerk for Trajectory Planning and Control", Robotica, Vol. 12, No.2, pp. 109-113.  [46]  Lai, J.-Y., Chen, H.-S., Lai, W.-B., 1995, "NC Cutting-path Generation From LaserScan Measurement Data", Journal of the Chinese Society of Mechanical Engineers, Vol. 16, No. 3, pp. 245-255.  [47]  Langeron, J. M., Due, E., Lartigue, C , Bourdet, P., 2004, "A New Format for 5-Axis Toolpath Computation Using B-spline Curves", CAD Computer Aided Design, Vol. 26, No. 12, pp. 1219-1229.  [48]  Lartigue, C , Due, E., Affouard, A., 2003, "Tool Path Deformation In 5-Axis Flank Milling Using Envelope Surface", Vol. 35, No.4, pp.375-382.  [49]  Lartigue, C , Thiebaut, F., Maekawa, T., 2001, "CNC Tool Path In Terms of B-Spline Curves"", Computer Aided Design, Vol. 33, No. 4, pp. 307-319.  [50]  Lin, R.-S., 2000, "Real-Time Surface Interpolator For 3-D Parametric Surface Machining On 3-Axis Machine Tools", International Journal of Machine Tools and Manufacture, Vol. 40, No. 10, pp. 1513-1526.  [51]  Lorenzo, S., Bruno, S., 1996, Modeling and Control of Robot Manipulators, McGrawHill Companies, New York.  [52]  Macfarlane, S., Croft, E. A., 2003, "Jerk-Bounded Manipulator Trajectory Planning: Design For Real-Time Applications", Transactions of IEEE on Robotics and Automation, Vol. 19, No.l, pp.42-52. .  [53]  Makino, H., Ohde, T., 1991, "Motion Control of the Direct Drive Actuators", Annals of CIRP, Vol. 40, pp. 375-378.  [54]  Matlab Reference Guide, 1992, The Math Works Inc.  [55]  Muller, M., Erdos, G., Xirouchakis, P., 2004, "High Accuracy Spline Interpolation for 5-Axis Machining", Computer Aided Design, Vol. 36, No. 13, pp. 1379-1393.  [56]  Park, H., Kim, K., Lee, S.-C, 2000, "Method for Approximate NURBS Curve Compatibility Based On Multiple Curve Refitting", Computer Aided Design, Vol. 32, No. 4, pp.237-252.  Bibliography [57]  115  Piazzi, A., Visioli, A., 2000, "Global Minimum-Jerk Trajectory Planning of Robot Manipulators", IEEE Transactions on Industrial Electronics, Vol. 47, No. 1, pp. 140149.  [58]  Piegl, L., Tiller, W., 1997, The Nurbs Book, Springer Verlag.  [59]  Pritschow, G., 1997, Course Notes: Steuerungtechnik der Werkzeugmaschinen und Industrieroboter (Control Techniques of Machine Tools and Industrial Robots), Institute of Control Technology for Machine Tools and Manufacturing Units, Stuttgart University, Germany.  [60]  Pritshow,. G., Rogers, G., 2004, "Decentralized CNC-Interpolation Using Intelligent Drives", University of Stuttgart, Institute for Control Engineering of Machine Tools and Manufacturing Units (ISW).  [61]  Qiu, H., Cheng, K., Li, Y., 1997, "Optimal Circular Arc Interpolation for NC Tool Path Generation in Curve Contour Manufacturing", Computer Aided Design, Vol. 29, No. 11, pp.751-760.  [62]  Rao, N., Ismail, F., Bedi, S., 1997, "Tool Path Planning for Five-Axis Machining Using The Principal Axis Method", International Journal of Machine Tools and Manufacture, Vol. 37, No. 7, pp.1025-1040.  [63]  Rehsteiner, F., Renker, H. J., 1993, "Collision-Free Five-Axis Milling of Twisted Ruled Surfaces", Annals of CIRP, Vol. 42, No.l, pp457-461.  [64]  Ruegg, A., Gygax, P., 1992, "Generalized Kinematics Model for three to Five-Axis Milling Machines and their Implementation in a CNC", Annals of CIRP, Vol. 41, No. l,pp.547-550.  [65]  Sarkar, B., Menq, C.-H.,  1991,  "Smooth-Surface  Approximation and Reverse  Engineering", Computer Aided Design, Vol. 23, No. 9, pp.623-628. [66]  Sarkar, B., Menq, C.-H., 1991, "Parameter Optimization In Approximating Curves and Surfaces to Measurement Data", Computer Aided Geometric Design, Vol. 8, No.4, pp.267-290.  [67]  Shpitalni, M . , Koren, Y., Lo, C C . , 1994, "Real-Time Curve Interpolators", CAD Computer Aided Design, Vol. 26, No. 11, pp. 832-838.  [68]  Simon, D., Isik, C , 1993, "Suboptimal Robot Joint Interpolation Within UserSpecified Knot Tolerances", Journal of Robotic Systems, Vol. 10, No.7, pp.889-911.  Bibliography [69]  116  Siemens A G , 2004, "Sinumerik 840D/840Di/810D Werkzeug und Formenbau Handbuch", Erlangen, Germany.  [70]  Siemens  AG,  2002,  "Advanced  Programming  Guide  SINUMERIK  840D/840Di/810D", Erlangen, Germany. [71]  Srinivasan, L.N., Ge, Q.J., 1998, " Fine Tuning of Rational B-Spline Motions", ASME Journal of Mechanical Design, Vol. 120, No. 1, pp.46-51.  [72]  Srinivasan, L.N., Ge, Q.J., 1996, "Parametric Continuous and Smooth Motion Interpolation", ASME Journal of Mechanical Design, Vol. 118, No. 4, pp.494-498  [73]  Trepanier, J.-Y., Lepine, J.L., Pepin, F., 2000, " Optimized Geometric Representation for Wing Profiles Using NURBS", Journal of Canadian Aeronautics and Space, Vol. 46, No.l,pp.l2-19.  [74]  Tsai, Y.-F., Farouki, R.T., Feldman, B., 2001, "Performance Analysis of CNC Interpolators for Time-Dependent Feedrates Along PH Curves", Computer Aided Geometric Design, Vol. 18, No. 3, pp.245-265.  [75]  Vassilev, T., I., 1996, "Fair Interpolation and Approximation of B-splines by Energy Minimization and Points Insertion", Computer Aided Design, Vol. 28, No. 9, pp. 753760.  [76]  Visioli, A., 2000, "Trajectory Planning of Robot Manipulators By Using Algebraic and Trigonometric Splines", Robotica, Vol. 18, No. 6, pp. 611-631.  [77]  Voliotis, S.D., 1992, "Orientation Planning in Continuous Path Applications for Wrist Partitioned Manipulators", IEE Proceedings, Part D: Control Theory and Applications, Vol. 139, No. 6, pp.495-502.  [78]  Wang, F - C , Yang, D.C.H., 1993, "Nearly Arc-Length Parameterized Quintic-Spline Interpolation for Precision Machining", Computer Aided Design, Vol. 25, No. 5, pp. 281-288.  [79]  Wang, F . - C , Schofield,  S., Wright, P.K., 1996,  "Real Time Quintic Spline  Interpolator For An Open Architecture Machine Tool", ASME Dynamic Systems and Control Division, Proceedings of the  1996  ASME  Engineering Congress and Exposition, Vol. 58, pp. 291-297.  International  Mechanical  Bibliography [80]  117  Wang, F . - C , Wright, P.K., 1998, "Open Architecture Controllers for Machine Tools, Part 2: A Real Time Quintic Spline Interpolator", ASME Journal of Manufacturing Science and Engineering, Vol. 120, No. 2, pp. 425-432.  [81]  Week, M . , Meylahn, A., Hardebush, C , 1999, "Innovative Algorithms for SplineBased CNC Controller", Production Engineering Research and Development in Germany: Annals of the German Academic Society for Production Engineering, Vol. 6, No. l,pp.83-86.  [82]  Wu, C. H., Jou, C C , 1989, "Planning and Control of Robot Orientational Path", IEEE Transactions on Systems, Manufacturing and Cybernetics, Vol. 19, No. 5, pp.12341241.  [83]  Xu, H.-Y., Dai, J., Tarn, H.-Y., Zhou, Y., 2003, "Angular Feedrate Interpolation for Three-Dimensional Implicit Curves", International Journal of Production Research, Vol. 41,No.l5,pp.3461-3478.  [84]  Yang, D.C.H., Golub, A.D., 1992, "Improved Trajectory Generation for Precision CNC Milling of Arbitrary Contours and Surfaces", ASME Production Engineering Division, Winter Annual Meeting, Vol. 52, pp. 105 - 117.  [85]  Yang, D.C.H., Golub, Alex D., 1994,  "Precision Trajectory Generation For CNC  Milling of Arbitrary Contours and Surfaces", International Journal of Machine Tools & Manufacture, Vol. 34, No. 7, pp. 1005-1018. [86]  Yang, D.C.H., Kong, T., 1994, "Parametric Interpolator Versus Linear Interpolator For Precision CNC Machining", Vol.26, No. 3, pp.225-234.  [87]  Yang, M.-Y., Park, J.-H., 2002, "A Study On An Open Architecture Cnc System With A Nurbs Interpolator For WEDM", International Journal of Advanced Manufacturing Technology, Vol. 19, No. 9, pp. 664-6683.  [88]  Yau, H.-T., Chen, J.-S., 1997, "Reverse Engineering of Complex Geometry Using Rational B-Splines", International Journal of Advanced Manufacturing Technology, Vol. 13, No. 8, pp. 548-555.  [89]  Yau, H.-T., Kuo, M.-J., 2001, "NURBS Machining and Feed Rate Adjustment For High-Speed Cutting of Complex Sculptured Surfaces", International Journal of Production Research, Vol. 39, No. 1, pp. 21-41.  Bibliography [90]  118  Yeh, S.-S., Hsu, P.-L., 2002, "Adaptive-Feedrate Interpolation For Parametric Curves With A Confined Chord Error", Computer Aided Design, Vol. 34, No. 3, pp. 229-237  [91]  Yeung, C.-H., 2004, M.A.Sc. Thesis: "A Three - Axis Virtual Computer NumericalControlled (CNC) System", University of British Columbia, Department of Mechanical Engineering, BC, Canada.  [92]  Zhang, Q.G., Greenway, R.B., 1998, "Development and Implementation of A NURBS Curve Motion Interpolator", Robotics and Computer-Integrated Manufacturing, Vol. 14, No.l, pp. 27-36.  [93]  Zhiming, X., Jincheng, C , Zhengjin, F., 2002, "Performance Evaluation of a RealTime Interpolation Algorithm for NURBS Curves", International Journal of Advanced Manufacturing Technology, Vol. 20, No.4, pp.270-276.  119  Appendix A Formulation of the Minimum Jerk Feed Profile The displacement profile for the k knot span segment was expressed in Chapter 4 as:  s (r)= A T +B T +C T +D T +E T 5  k  4  k  3  k  +F  2  k  k  k  ,0<r<T  k  (A.l)  k  In order to avoid any ill conditioning, the time parameter z and segment durations can be normalized by the longest expected segment duration T  max  at the beginning, resulting in the  following feed spline equation,  s (v) = a v +b v 5  k  where  +c v  4  k  v = r/T  +d v  3  k  2  k  k  k  ,X = T /T  max  k  k  +e v + f  and  max  , 0<v<^  k  a = T* A ,b k  ax  k  k  (A.2)  = T* B ,...,f ax  k  k  =F . k  The  normalized velocity, acceleration, and jerk equations of the given quintic time spline are as follows:  s (v) = (1 /  ) [5a v + 4b v + 3c v +2d v + e ] 4  k  3  k  2  k  k  k  k  ^(v) = ( l / 0 [ 2 ( H v +12V +6c v + 2d ] 2  3  k  s {v)=(l/Ti )[60a v  2  k  x  k  . (A.3)  k  +24b v + 6c ] k  k  The local minimum jerk objective function was given in Eq. (4.35) (Chapter 4) is expressed using the normalized time spline parameters as:  3600a v + 2880a,Z>,v + 2  J  Vk  T  5  J [(l20a c k  4  3  + 576^ )v + 2SSb c v + 36c 2  k  \dv.  2  k  k  The integration in Eq.(A.4) results in the following matrix formulation:  (A.4)  Appendix A. Formulation of the Minimum Jerk Feed Profile  1204 0 72% 0 1204 72^ 36\ 0 0 0 0 0 0 0 0 0 0 0 0 0  where 9 = [a k  b  k  f] k  dropping the l/T^  36<HJ  3604  192^  0 0 0 0 0 0  0 0 0 0 0 0  a  k  (A.5)  »  fk v  •  is the normalized parameters of the displacement spline. By  T  k  7204  120  at the beginning of Eq.(A.5) and replacing with 1/2 for mathematical  convenience, the minimum jerk objective for the k feed segment can be rewritten as: th  Vk  '  J  K  k { \ ) -  f  t  k  ~ ~ j f k  '  K  k  (A.6)  '&k  The overall minimum jerk objective in Eq.(A.6) can then be stacked for all the segments, where the local parameters and weighting matrices for local segments are 6 and K as: k  1  Ns-l  R  k=Q  T  0!  Jv=HJk=-\0  1  k  0o 0,  Ns-\ \'  9Ns-l. *6Nsxl  Position (C°), C , C and C boundary conditions are imposed at the segment junctions to 1  2  3  connect the time spline segments together over the complete toolpath. For the k segment, th  the following position boundary conditions must hold,  (A.8)  Appendix A. Formulation of the Minimum Jerk Feed Profile  121  where S is the arc length of the k' segment. The boundary condition is rewritten in the k  matrix form as:  0  0  0  0  o  Eft  1 I  (A.9)  The knot position constraints from all segments can be grouped together as:  L°  ' t°o ' (A.10) 0Ns-l L"(i.)  Following that, the jerk continuous feed profile is generated by imposing the velocity, acceleration, and jerk continuity constraints. For instance, at the junction of k  th  and k+l  st  segments the boundary conditions are given in Eq.(4.32) (Chapter 4) as,  A+1  3  (0)1  ^(A)-^ .(o)=o  (0)  •s (^)-S ,(0) = 0 t  k  k+  where k = 0,l,...,N  s  (A.11)  t+  s {\)-s {0)  s {\)=s M, k  +  M  = 0\  - 2. Substituting the expressions for for s ,  (A.3) into Eq. (A.l 1), the matrix equation is written,  k  s  k  and s' , from Eq. k  Appendix A. Formulation of the Minimum Jerk Feed Profile  122  44  34  24  1 0  0 0  0  0  -1  0  204  124  64  2  0 0  0 0  0  -2  0  0  = 0  604  244  6  0  0 0  0 0 -6  0  0  0  0  5K  "0" (A.12)  4 (4) 23  The derivative boundary conditions for N -\  segment junctions (Eq. (A.12)) are  s  grouped together and stacked in the following matrix equation as:  ~L  0  123  •^0  U  0 U  3x6  L  3x6  0,  1 2 3  (A.13)  ^1 123 Afc-2  0  0,  r  J  (i,„„), acceleration (s  If initial conditions of velocity  ),  and jerk (s' ) are provided,  injt  they can be included into the formulation as, s(o) 0  in!t  = s , s(o) = s,.„,., and '^(0) = 's;m7  0  mil  in the  feed profile. The boundary conditions expressed using the Eq. (A.3) are given in the following matrix equation as:  "0 0 0 0 1  0,  0 :  o  0 0 0 2 0 0 ; 0  0 6 0 0 0 ;  0, w  N  T s' T s 3  0,  Similarly, if final kinematic conditions of velocity (s ), fiml  ('sfmai)  a r e  (A. 14)  2  )  acceleration (s' ), and jerk final  provided at the end of the toolpath, they can be included as s _ (A, _ ) = Ns  final  S  '  Ns-l  S  fmai  s  a  n  d  written in the matrix form:  {  Ns  }  s , final  Appendix A. Formulation of the Minimum Jerk Feed Profile  5  03x6(/Vi-l)  4,-i  34,-i  44-i  204,-1 124,-1 604-i 244-,  6  O  0 0 0 0 0 1  24,-i 2  4,-i 6  123  0  0  0  ^max^f  T s i ... 2  (A.15)  T  il {^Ns )  At the segment junctions; position, I , 2 , and 3 derivative continuity constraints are s  n  r  given in Eq.(A.9) and Eq.(A.12). The initial and final boundary conditions are provided for the feed profile from Eq. (A. 14) and Eq. (A.15). All of the constraints are stacked together in the following matrix equation form:  / (4 l23  L final (4,-1  •6  L(X)-0  (A.16)  =£  )  The minimum jerk optimization problem can now be expressed using Eq. (A.7) and Eq. (A.16) as:  min -6 -K(X)-0  subject to: L{X)-0  T  The  (A.17)  above constrained quadratic minimization problem can then be solved using the  Lagrange A = [AQ  =£  2  0  Multipliers A  {  ...  A ], T  5NS+2  method.  Introducing  vector  of  Lagrange  multipliers  the objective function can be expressed as,  j (0, A)=- 0 • K{X) -0+A T  v  the  T  -0-4]  (A.18)  Appendix A. Formulation of the Minimum Jerk Feed Profde  By equating the partial derivatives to zero ( /QQ V  =  u  a n  124  d  /QA  V  =  ^ y^ds the linear  equation system:  K{X)-0  + L {A)-A  L{X)-0  =Z  T  = 6\  ~K(X)  J  L (X)  ~e  0  A  T  (A.19)  which has full rank as long as all the segment durations are nonzero. The solution of the above matrix equation yields the coefficients (0) of the minimum jerk feed profile. The coefficients are de-normalized as A = 1 / 7 ^ -a ,B k  profile is generated using Eq.(A.3).  k  k  = \/T*  max  •b ,...,F k  k  =f  k  and the feed  125  Appendix B Computation of the Gradients in Feedrate Optimization B . l . Gradient of the Objective The minimum time objective was given in the Eq. (4.36) in Chapter 4 as:  Ns-l  k=0  and in normalized form of  ^ = Z^=^Z4 Ns-l  Ns-l  k=0  k=0  (B.l)  Hence, its gradient with respect to the normalized time duration of the k segment in the th  B-spline representation X is expressed as: k  ^r =T  dx  T  (B.2)  max  y  }  The gradient with respect to the vector of all normalized segment durations X = [X  0  X  x  ...  X _]  ff  T  Ns  x  is obtained in a vector form as:  =^ - [ l  1  l l *  (B-3)  Appendix B. Computations of Gradients in Feedrate Optimization  126  B.2. Gradient of the Feed Profile Coefficients The vector of the Quintic displacement profile coefficients: 0 is obtained by solving Eq. (A. 19) in Appendix A as:  K(X) • 0 + L {X) -A = Q\ ~K(X) T  L {X)  _L(X)  T  o  ~e  "o"  A  A.  (B.4)  It is known that the b vector on the right side of Eq. (B.4) is a constant depending only on the toolpath arc lengths and the given boundary conditions. The differentiation of both sides of the Eq.(B.4) with respect to the design variables /I • yields the following equation:  dA cU.  X  (B.5)  A and X are obtained from the solution of Eq.(B.4), the gradient expression dA/dXj is constructed by considering the definitions of the K(X)  and L{X)  matrices given in Eq. (A.7)  and (A. 19) (Appendix A). Basically, the expression for dX/BXj is obtained by solving the linear equation system in Eq.(B.5), and the gradient of the feed profile coefficients dd/dXj is extracted from:  dX dX,  = [de /dXj T  dA /dXjj T  (B.6)  The gradient matrix, consisting of partial derivatives of the feed profile coefficients with respect to all normalized segment durations is constructed as:  Appendix B. Computations of Gradients in Feedrate Optimization da /dZ  0  da 1cU,  dbJdX  Q  db /dA  dfo/ 4 da /dA  .  0  df/dA, daj8A  db /dA  0  dbJdA,  . ••  8fi/dA  df/dA,  . ••  0  0  0  8  x  x  36  127  ••  l  0  ••  l  0  b /dA _  d  Na  aJdA .  d  Na  d°JdA _ Na  x  fl ^Na-l  d  d  Ns-J K-i  da  Sb .JdA  db _JdX  f s-\l 4  df -Jd^  Ns  0  Ns  d  d  N  d  db _JdA _ Ns  x  Ns  {  f s-J ^Ns-\  d  d  N  Ns  16NsxNs  (BJ)  B.3. Gradient of the Kinematic Profiles The arc displacement, feedrate, tangential acceleration, and jerk profiles in the k  tF  segment have been given as:  s M = a v + b v' + c v + d v +e v + f 5  k  3  k  k  2  k  k  ^(v) = ( l / ^ ) ( 5 f l v + 4 V 4  ^(v) = ( l / ^ ) ( 2 0 a , v s {v)=(l/TL)(60a v  2  k  k  k  2  t  3  k  +3c,v + 2d v + e )  3  k  + I2b v +6c v 2  k  k  k  + 2d ) k  (B.8)  +24b v + 6c ) k  k  where 0 < v < X . The normalized time parameter v is redefined as: k  v = A p where 0 < p < 1 k  (B.9)  Appendix B. Computations of Gradients in Feedrate Optimization  128  The displacement profile is then written as:  h  = "XM  + bXf* + XS  5  + d£p  C  (B.10)  + eXp + f  2  k  k  k  k  k  The gradient of the displacement profile in the k segment with respect to the vector of normalized segment durations X = [X  /\  Q  ...  /\ _ ] is obtained as: Ns  x  (B.H)  where the gradients da /dA ,  db /dA ,  T  T  k  gradient of the k durations is  th  df /dX  T  k  are obtained form Eq.(B.7). The  normalized segment duration with respect to all normalized segment  dA /d/l  T  k  k  = [0  ... 0 1 0  ... 0]  where only the k  th  entry is one.  Considering the definition of s in terms of normalized parameter in Eq.(B.8) and factoring out the p term out of the expression in Eq.(B.ll), the gradient of the displacement profile with respect to the normalized segment duration is obtained as:  f + ^F +  ^T J {v) m  k  dX  k  (B.12)  Following a similar procedure, the gradients of the feedrate, tangential acceleration, and jerk profiles are also obtained as:  Appendix B. Computations of Gradients in Feedrate Optimization da  5v  dX  +4v  kk  4  T  3  dX  T  db  3v  k  2  — f r +  dX  T  dc  k  - r ^T +  dX  dd  de  k  2v  dX  54  2  k  max k  dX  T  "  T  129  S  5^  r  dc da. , - db ,„ dc „dd _ ... / N5/L 20v —%• +12v — + 6v — - + 2 — + uT's I v)— dX dX dX QJT A1 t~ max k\ J dX 3J  2l  kk  T  kk  T  k  3  k  T  k  T  J  k  T  (B.13)  T  6 0 v % + 24v %  6 %  2  dX  T  +  +  <  ^(v)%  B  B . 4 . Gradients of the Constraints The velocity, acceleration and jerk constraints have been given in Chapter 4, Eq.(4.60) as:  -s(t)lf -\<Q,  C- (t) =  mca  f  = X(t)/v  -\<0,C (t)  Xmax  Ymax  cKo-= Z(t)lv  Zmm  max  = -X(t)/v  x  c;ro== Y(t)/v  -s(t)lf <0 -1 <0  Xmax  -1<0,  C (t) = -Y(t)/v  -1 < 0,  CTfO - - Z(t)/v  Y  -1<0  Ymai  = 0,f 0 / ^ _  -1<o,  c-jt) = -e (t)/v  = Oc(0/v  -1 < 0,  C: (t) = - 9 (t)/v  6c  max  C (t) =x(t)ja  -1 < o,  = Y(t)/a C y(0-  -1 < o,  C\(t)- = Z(t)/a  -1<0,  x  Xmax  +  Ymax  Zma  A  c  0c  = - wK m a t -iso  C~(t) = -Z(t)/a  - -^  c  1  -1 < o,  c~ jt)=- e (t)/v  C (t) = X(t)lj  -1 <0,  C' (t) = -X(t)/j  — 1 < 0,  Cy(t) = - Y(t)lj  -1<o,  c=ro=-z'ro/«z  c  ec max  Xmax  x  rO Y (t)jj =  Ymax  c|ro=z(t)/a  Zmax  x  Bc  -1 <o -1 < 0  Ymax  mai  C~ (t) = - e (t)/ ,  c max  -1 < o,  q t)=- e (t)/v,  ci (t) = e (t)/v,  Upper Bound Constraints  9a  c(  A  c  0  - I<o  max  Xmax  -1 < o,  A  c  c  max  ci 0) = e (ofa A  0  Acceleration,  -1<0  Zmax  c; (t) = e (t)/v c  -1 < 0  max  Q C ; = - i r o A ^ -iso  =<?,r o K - - * °- ^ ^=- ^ 1  -I<o  6Amax  C  Velocity,  -1 < 0  Zmax  v A  -1so max  c max  Lower Bound Constraints  -1 < o -1 < o  Jerk Constraints  (B.14)  Appendix B. Computations of Gradients in Feedrate Optimization  The feedrate upper bound constraints C {t) = s(t)/f +  f  max  130  -1 < 0, its gradient with  respect to normalized segment durations s obtained as:  OC}  1 ds{t)  =  ^  (B.15)  La*  where the gradient of the feedrate profile  is given in Eq.(B.13). The axis kinematic  profiles have been given in Chapter 4, Eq. (5.58) as:  M(t)  =  M(t)  =  M(t)  = M -sitf  where  M -s{t) s  M .s(t) +M -s(t) 2  ss  (B.16)  s  +3M -s(t)-s(t)+M -sit)  sss  SS  M(s) = [x(s)  M  Y(S) Z(S) 0 (S) a  0 (s)J.  The X axis velocity upper bound  c  constraint can be obtained by combining Eq.(B.14) and (B.16) as:  Q = —  * « - !=—  X (s)-s(t)-l<0 s  (B.17)  and its gradient with respect to normalized segment durations is obtained as:  dC _ x  1  ./ \ ds  ds  (B.l 8)  The X axis acceleration upper bound constraint is obtained from Eq.(B.14) and Eq. (B.16) as:  C\ = —  Xit)-1 = —  [x {s). s{t) + X, • s{t)]-1 < 0 2  ss  (B.19)  Appendix B. Computations of Gradients in Feedrate Optimization  131  and its gradient with respect to normalized segment duration is obtained as:  dl  (X j(tf  + Xj(t))^  s  a.  T  +  2Xj(t)^r  X  +  s  §  (B.20)  The X axis jerk upper bound constraint can also be obtained from Eq.(B.14) and Eq. (B.16) as:  C =• +  1  ^t)-l  = ^[xJs)-s{ty+3X -s{ty-4t)+X -s{t)}-\<0 ss  X  (B.21)  s  A  and its gradient with respect to the normalized segment durations is expressed as:  dC _ x  1  (* J  3  SS  + 3X Js s  + X s)-^r ss  +  3(X J> S  +X  „\ a?+.,3XJ „ .| aL? +. „ X. ds  j ) - ^ r  (B.22)  The gradients of the lower bound constraints have identical expressions to those of the upper bound constraints, with opposite sign (Eq.(B.14)). The gradients for the Y, Z, A and C axis constraints can be obtained by replacing the occurrences of 'X' with the rest of axis variables in Eq. (B.l8), (B.20), and (B.22).  

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0080755/manifest

Comment

Related Items