Numerical implementation of drop spin and tilt method for five-axis tool positioning for tensor product surfaces

This paper presents an extension of multipoint machining technique, called the Drop Spin and Tilt (DST) method, that spins the tool on two axes, allowing for the generation of multiple contact points at varying distances around the first point of contact. The multiple DST second points of contact were used to manually generate a toolpath with uniform spacing between the two points of contact. The original DST method used a symbolic algebra package to position the tool on a bi-quadratic surface; our extension is a numerical solution that allows positioning a toroidal tool on a tensor product Bézier surface, which we tested on bi-degrees 2, 3, 4, and 5. On the three bicubic surfaces we tested, the average number of initial seeds for our numerical method was less than 5 across all three surfaces, and the average number of Newton iterations was less than 20 except when the curvature of the tool closely matches the local surface curvature (in which case more iterations were required). Furthermore, we investigate the spread of possible second points of contact as the tool is spun around these two axes, demonstrating the feasibility of using the method to control the machining strip width.


Introduction
5-axis machining is a reliable and efficient method to machine geometrically complex parts found in tool and die manufacturing, aerospace, automotive, and marine industries. 5-axis machines may use ball, flat, and toroidal (radiused) end milling cutters, but the cutting characteristics of the toroidal cutters are better at the bottom tip of tool [1] Sandeep Kumar  University of Waterloo, Waterloo, ON, Canada thus offering better surface finish. The geometrical details of a toroidal end milling tool are shown in Fig. 1.
5-axis NC machining offers improved cutting efficiency, such as higher material removal rate and surface finish, over 3-axis machining. These advantages are drawn from the flexibility in the machine kinematics that enables tool positioning in a close proximity to a curved design surface.
Irrespective of the growing complexity of the parts, the focus of most current research has been on improving the machining times while maintaining the surface finish. Several 5-axis tool positioning strategies for machining of sculptured surfaces have been developed by researchers to address these issues. These approaches can be classified into two sets.
In the first set, the tool is positioned so that tool curvature matches the local curvature of the design surface in the area surrounding the single contact point. The curvature matching methods do not guarantee a gouge free toolpath as such. Thus, gouge detection is integrated with the positioning and inclination as a secondary step.
In the second set of solutions, the tool simultaneously touches the workpiece at more than one contact points and thus is termed as Multipoint Machining methods (MPM) [2]. This multipoint contact is accomplished by using enhanced tool positioning and orientation algorithms [3]. Machining times can be significantly reduced by having such tool positions and orientations that yield wider machining strips while meeting the surface finish constraints. In such methods, the first point of contact is relatively easy to obtain, while the second point of contact is more challenging. Several methods have been suggested by researchers for finding the two points of contact; Liang et al. [4] presented a survey of 5-axis tool positioning methods that included an overview of such techniques.
The first MPM technique of interest is the Drop and Tilt Method (DTM), which is a two-step procedure to find gouge-free, 5-axis tool positions on curved surfaces. Duvedi et al. [5] developed DTM for triangulated curved surfaces and later extended it to Bézier tensor product surfaces [6]. Although DTM is computationally robust, it cannot control the distance between the two contact points for a given location on a curved surface. Sharma et al. [7] took DTM further and developed the Drop Spin and Tilt Method (DSTM), which is a three-step method that involves the merger of DTM and DSM [8]. The DSTM method was developed mainly to regulate the distance between the two contact points for given spin angle to generate the 5-axis tool position. Thus, DSTM addressed the key limitation of DTM and enabled the control over the machining strip width by incorporating an intermediate tool-spin step.
As a proof of concept, the DSTM was tested on bi-quadratic Bézier tensor product surfaces within the Maple TM symbolic algebra environment. However, even with bi-quadratic test surfaces, the solution involved highorder transcendental rational equations. The solution of these equations was reasonably challenging in the symbolic algebra environment. Therefore, dealing with higher order surfaces using symbolic solver tools is expected to be even more computationally expensive and infeasible for typical computing facilities available in industry.
The focus of this work is to develop an algorithm to make DSTM a system independent, stand alone application that can be reliably used to generate gouge-free 5-axis tool positioning of radiused end mill tool for machining of higher order tensor product surfaces. The proposed DSTM numerical algorithm has been implemented in C++, and tested on bi-quadratic, bi-cubic, bi-quartic, and bi-quintinc Bézier surfaces. In this paper, we give details on the results for the bi-cubic surfaces to allow for comparison to earlier work. We also perform a set of tests to show the feasiblity of using the method to control the machining strip width.

Background
Most industrial products like dies and molds used in manufacturing of sheet metal parts, automotive and space applications, have curved surfaces. Such curved shapes are generally modelled by Bézier or B-spline surface patches. In the present work, MPM toolpaths for machining a tensor product Bézier patch have been developed by numerical implementation of DSTM.
For generating the toolpath, the Bézier patch is projected on the z = 0 plane, also called the tool drop footprint plane. A zig-zag pattern describing the motion of the tool is developed within the rectangular region. Typically, the zig-zag pattern is discretized using two parameters, namely the feed forward step (f-step) and the side step (s-step), as shown in Fig. 2. The discretized points are termed as tool drop positions. The toolpath is generated by moving the tool between such consecutive tool drop positions (T) and applying the three-step DSTM procedure to determine gouge-free tool positions for each position T. The machine controller drives the tool between the sequence of tool positions in a toolpath.
To find a single tool position, the DSTM method uses three steps: tool-drop, tool-spin, and tool-tilt. Note that during machining, the machine does not perform these In the first step, a radiused end mill cutter is vertically dropped on a surface S with the tool axist passing through a given tool drop position T as shown in Fig. 3a. The dropped tool finds the gouge-free highest point of contact on the surface S. This is the first point of contact P as shown in Fig. 3a. The point P lies on a circular pseudo-insert that represents the cross-section of the circular cutting edge of the toroidal tool. The surface normal through P isn p .
In this configuration, the tool can be rotated about the axis of the pseudo-insert, or about the surface normaln p at the contact point, without losing tangential contact with the surface. Multi-point methods have typically focused on rotating the tool about the pseudo-insert axis; however, in DSTM, the dropped tool is first rotated about the surface normaln p and then rotated about the pseudo-insert axis to get the second point of contact. This additional rotation about the surface normaln p provides a variety of possible second point of contact lying at a desired distance from the first point of contact.
In the second step of DSTM, the dropped tool is rotated about the surface normaln p by an angle α, while remaining tangent to the surface S at point P (Fig. 3b). The tooldrop axist generates a conical surface of revolution when the tool is rotated about the spin-axisn p . The spin-angle α determines the gap between the two points of contact, thus eventually controlling the cut strip width. A suitable value of spin-angle is chosen, which in turn is used to find the orientation of tool-axist in the tool-spin step.
Finally in the third step of DSTM, the tool is rotated about the tilt-axisî 3 , which passes through the center of the pseudo-insert and is perpendicular to the plane of pseudoinsert as shown in Fig. 3c. The tool is tilted towards the surface until the tool touches the surface at a second, gougefree tangency point, Q. This gives a gouge-free multipoint tool position so that tool is tangent to the surface S at two distinct locations, P and Q. The gap between the two contact points is controlled by varying the spin-angle α in the spin-step.

Mathematical model of drop spin and tilt method
The surfaces used in this work are tensor product Bézier surfaces. The representation of a Bézier tensor-product surface is given by with 0 ≤ u, v ≤ 1, where P i,j are control points, and B n i (u) and B m j (v) are the Bernstein basis functions given by We assume the xy-coordinates of the control points of the Bézier surface to be uniformly distributed in the rectangular footprint domain. This restriction conveniently represents the Bézier surface S(u, v) in the representation z = f (x, y). This facilitates the evaluation of the z-coordinate of the surface S, which simplifies the gouge checking procedure.
In this section, we derive the mathematical details of our method; additional implementation details along with pseudo-code are given in Appendix A.

Step 1: Tool-drop
In the first step, a toroidal tool is dropped vertically on the Bézier surface S(u, v), along tool-axis direction t = [0, 0, 1] while tool-axis also passes through footprint location T as shown in Fig. 3a. The toroidal tool is centered at T c , with radius of insert R i and offset radius R o . The location of the center T c of the dropped tool is given by where h represents the height of the tool center T c above footprint point T. The vertically dropping toroidal tool touches the surface S at point P = S(u 1 , v 1 ) as shown in Fig. 3a. From the Drop and Tilt Method [9], we know that the first point of contact S(u 1 , v 1 ) for a toroidal tool dropping alongt through footprint point T on surface S is wheren p andn r are the surface normal at point P and a projection of the surface normal perpendicular to the tool axist: Equation 2 can be solved for (u 1 , v 1 , h) to determine the first point of contact P , surface normaln p andn r . Alternatively, given u, v, we can determine the contact point P , surface normaln p andn r , an approach that we used in our implementation (Appendix A). Furthermore, the location of the center of the pseudo-insert O 1 through the first point P is A coordinate frame F 0 = {î 0 ,ĵ 0 ,k 0 } at O 1 is used to define the initial position of the dropped toroidal tool centered at T c . This frame is F 0 = î 0 ,ĵ 0 ,k 0 = n r ×t n r ×t ,n r ,t if the first point of contact P is on the outside of the toroidal tool; if P is on the inside of the torus, a similar frame is used but withî 0 =t ×n r . The dropped toroidal tool, T tor0 , with its center located at T c , is defined in frame {î 0 ,ĵ 0 ,k 0 } as, where θ and φ are the angles used to model the toroidal tool in the frame {î 0 ,ĵ 0 ,k 0 }. The orientation of the direction vectors for frame {î 0 ,ĵ 0 ,k 0 } is shown in Fig. 3a.
The drop method has been implemented and tested both in a symbolic algebra environment and in C++ [9].

Step 2: Tool-spin
In the second step, an orthogonal coordinate frame F 1 = {î 1 ,ĵ 1 ,k 1 }, called the tool-spin frame, is defined at O 1 as shown in Fig. 3b. The tool-spin frame is created by rotating the tool-drop frame F 0 aboutî 0 so thatk 1 is aligned with the surface normaln p , while the dropped tool stays fixed in the frame {î 0 ,ĵ 0 ,k 0 }.
The direction vectors of the tool-spin frame F 1 are The representation of the dropped toroidal tool with respect to frame F 1 can be obtained by using a rotation R n m , where R n m is a rotation matrix that maps a point from frame-n to frame-m: All subsequent rotations ensure that the tool definition frame F 1 is mobile and it rotates along with the tool in the tool-spin and tool-tilt steps.
In the tool-spin step, the dropped toroidal tool T tor1 , defined in frame F 1 , is rotated about the spin-axisk 1 passing through O 1 by spin angle α. The desired tool-axis orentation t is obtained by using a specified spin-angle α, which eventually determines the spacing between the points of tangency between the work surface and the toroidal tool in the third step.
The post-spin orientation of the toroidal tool is defined in frame F 2 = {î 2 ,ĵ 2 ,k 2 }, which is initially assumed to be coincident with F 1 and rotates along with the tool in the tool-spin step as shown in Fig. 3b. Frame F 2 is the tool-tilt frame, and its direction vectors are

Step 3: Tool-tilt
In the third step, the tool in the tool-tilt frame F 2 is further rotated by an angle β about the tilt-axisî 3 until it is tangent to the surface S at a second point of contact Q as shown in Fig. 3c. The resulting tool orientation is defined using another orthogonal coordinate frame F 3 = {î 3 ,ĵ 3 ,k 3 } with its origin at O 1 . The frame F 3 is initially coincident with frame F 2 , and rotates along with the tool through angle β in the tool-tilt step. The direction vectors of the F 3 are given as In the tool-spin and tool-tilt steps, the toroidal tool continues to maintain tangential contact with the surface S at P . After the spin and tilt rotations, the tool has a second point of contact with the surface S at point Q as shown in Fig. 3c.
The position and orientation of the rotated toroidal tool with respect to the standard Cartesian frame {x, y, z} given as T tor are Similarly, the orientation of the rotated tool-axist and the rotated major circle of torus C tor (centered at T C ) in frame {x, y, z} is given aŝ When the tilted toroidal tool touches the surface S at Q = S(u 2 , v 2 ), the pseudo-insert through Q has its center on the major circle of the rotated torus C tor . Thus, wheren q is the surface normal at S(u 2 , v 2 ) and is defined similar ton p . For a given spin-angle α, Eq. 6 is a simultaneous rational equation with three component equations in x-, y-, and zcoordinates but has four unknowns u 2 , v 2 , β, and φ. An additional constraint that ensures that the shortest distance between the surface normaln q through the second point of contact S(u 2 , v 2 ) and the rotated tool axist is zero is used to solve for the four unknowns: We solved Eqs. 6 and 7 for (u 2 , v 2 , β, φ), and in turn the values of these four parameters are used to compute the point Q, surface normaln q , and direction vectors {î 3 ,ĵ 3 ,k 3 }. Frames F 1 , F 2 , and F 3 , and the transformation matrix R 3 0 are used to determine the rotated toroidal tool T tor , tool axist , tool center T C , and second point of contact Q. The resulting tool positions from the numerical solver may gouge the surface. Thus, a gouge checking method [9] was implemented to ensure that the accepted tool positions do not gouge the surface beyond user-specified tolerance, while some gouging may be acceptable for actual machining, requiring no gouging seemed appropriate for our tests.
Pseudo-code for our implementation for this method are given in Appendix A.

Numerical implementation
The numerical DSTM algorithm computes the toolpath for multipoint tool positions for a tensor product, Bézier surface S (degree N × M) with uniformly spaced control points as shown in Fig. 4. Pseudo-code for our algorithm appears in Fig. 5. The toolpath is generated within the limits of X max and Y max along the X-and Y -directions, which also define the span of the control points of the uniformly spaced Bézier surface S. A zig-zag tool-drop trajectory is defined at z = 0 within the X max and Y max span, and is discretized into X Steps and Y Steps along the X-and Y -directions by using a feed forward step of f step and a side step of s step .
A toroidal tool, defined by a circular insert of radius R i placed at offset radius R o , is dropped vertically on such discretized locations (marked as T in Fig. 3a) to find a nongouging point of tangency with the surface. Subsequently, the dropped tool is turned about the surface normaln p Newton's method is used in the tool-drop and tool-tilt steps. To make the numerical solution independent of the degree of the Bèzier surface, the partial differentials of the surface and the surface normals are computed using the Richardson extrapolation. In the tool-drop step, Eq. 2 is solved with Newton's method to compute the first point of contact S(u 1 , v 1 ), and the height of the tool center for the dropped tool. A random initial seed is provided to the Newton's method to start the solution search. If Newton's method does not converge within a maximum stipulated number of trials, then another random initial seed is tried. If the tool-drop converges, then the solution is further tested for the gouge-free tangential contact between the tool and the work surface. If the converged solution yields a gouging tool position, then Newton iterations are aborted and another random seed point is reinitialized.
The dropped tool is then rotated about the surface normal n p by an spin angle α and is subsequently tilted about thê i 3 axis through an angle β so that the tool has a gougefree tangency at a second point of contact Q = S(u 2 , v 2 ). Newton's method is used to compute the tilt solution and the search starts with a random initial seed. In case the tooltilt solution for the previous/adjacent foot print location is available, then the previous solution is used as the initial seed; if Newton's method does not converge with this initial seed, then in the next trial, a random initial seed is used. If Newton's method converges, then the computed tooltilt solution is tested for its proximity with first point of contact, tool-axis inversion, and a gouge check. These three checks ensure that the first and second points of contact are separated by a user specified gap vicinity , the tilt angle β is negative, and the tilted tool is not gouging the surface beyond the permitted gouge tilt tolerance ( Gouge ). If the tool-tilt solution fails on any of these three checks, Newton's method is reset with a new random initial seed. In case no valid tool-tilt solution is found within a maximum number of initial seed solutions (N tilt ), then it is construed that for the given spin angle α, the tool tilt is infeasible and the first and second points of contact are taken to be identical. Details of the numerical implementation are available in the Appendix A.

Test surfaces
We implemented the numerical solver for DSTM algorithm in C++ and tested it on bi-quadrtic, bi-cubic, bi-quartic, and bi-quintic Bézier surface patches. Here, we give details on three different bi-cubic test surfaces, concave, convex, and  saddle; the control points used for these patches are given in Tables 1, 2, and 3. The x-and y-coordinates of the control points of all the test surfaces are distributed uniformly over a span of 150mm×150mm as shown in Fig. 4. We choose to give the results of these three surfaces in this paper for comparison to earlier work [9]. The tool-drop footprint surface was obtained by projecting the surface patch on the z = 0 plane. The zig-zag toolpath was generated on the tool-drop footprint plane. This toolpath is composed of 16 machining passes with a side step of 10 mm. In each machining pass, a total of 61 tool positions were computed with a feed forward step of 2.5 mm between two consecutive tool-drop positions. Therefore, in this DSTM study, the toolpath for each test surface was composed of 976 discretized tool positions. The toroidal tool with an offset radius R o =6.7 mm and insert radius,R i =6.0 mm was dropped vertically on each discretized location as discussed in Section A.2.
Newton's method was used to obtain the tool-drop and tool-tilt solutions. The partial derivatives required in the computation of the Jacobian matrix were computed using the Richardson extrapolation, with an initial step size (h step ) of 0.001 mm and a convergence tolerance ( Richard ) of 10 −5 mm. The seed points required to initialize Newton's method were randomly generated.
In the tool-drop step, tolerances of Drop = 10 −6 mm was used as convergence criteria in Newton's method. Each solution for first point of contact from tool-drop step was examined using gouging check algorithm with a gouge tolerance of Gouge = 10 −3 mm. The maximum number of random seed points (N Drop ) and the maximum Newton iterations (N I ter ) in the tool-drop step was kept as 500 and 150 respectively. Furthermore, the dropped tool is rotated by spin-angle (α) about surface normal (n p ) passing through first point of contact at P . In this work, the toolpaths were generated with values of spin-angle varying from 0 to 180 • with uniform increments of 15 • .
Finally in the tool-tilt step, the number of random trial seed points (N T ilt ) and number of Newton iterations (N I ter ) was taken as 50000 and 150 respectively. A relatively higher value of N T ilt was employed as a preventive measure so that a feasible tilt-solution is computed at each drop location. For the tool-tilt step, tolerances of T ilt = 10 −4 mm and Gouge = 10 −3 mm were used as convergence criteria in Newton's method and for gouge check respectively. The tolerances used for convergence in tool-drop and tool-tilt steps are chosen to ensure that the resulting 5-axis tool positions are accurate in the range of 10 −3 mm.
DSTM algorithm statistics are given in Table 4. This table shows the maximum and the average number of the random seed points, along with maximum and the average number of iterations taken for Newton's method in the tool-drop and tool-tilt steps for all three test surfaces corresponsing to various values of spin-angle.
The maximum number of randomly generated initial seeds taken for the tool-drop step across all surfaces, at different values of spin-angle α, is less than 100. Moreover, the average number of initial seeds is also observed to be less than 5 when averaged over 976 tool positions for each spin-angle. These results indicate that in a tool path (comprising 976 tool-drop locations) except for few hard tool-drop locations, the majority of the locations need only a few initial seeds to converge to the acceptable solution. Similarly, the maximum number of Newton's method iterations in the tool-drop step is found to be less than 20, except for the saddle surface. In the tool-drop step, the maximum value of N I ter for the saddle surface is relatively higher when compared to the other test surfaces. This may be attributed to the inflection regions on the saddle surface where the curvature transitions exist. However, the average number of Newton iterations for the saddle surface is marginally higher as compared to the other two surfaces.
The tool-tilt step of the DSTM algorithm is computationally more expensive than the tool-drop step. This extra expense is indicated by the relatively higher values of the    maximum number of initial seeds used in the Newton solver in the tool-tilt step. Specifically, in the tool-tilt step for the saddle surface, at spin-angles 15 • , 30 • , and 150 • , the maximum number of initial seeds used was exceptionally higher (greater than 1000). This need for more initial seeds is attributed to the conditional check incorporated in the tool-tilt step to determine the second point of contact Q away from the first point of contact P . This spacing between P and Q is ensured by specifying a minimum parametric spacing between the two contact points P and Q using the where vicinity = 0.001. However, when the curvature of the tilting tool closely matches the local surface curvature, there is hardly any scope for the tool to tilt because as the tool is tilted, the tool curvature decreases, and to have a gouge free tool position, the tool curvature can not be smaller than the local maximum surface curvature. In such situations, the second point of contact tends to lie close to the first point of contact and thus it is hard to have a large gap between the two points of contact. In such cases, when the second contact point is not achieved within first 1000 initial seeds, the tolerance for gap between two points of contact, given by vicinity , is relaxed to 0.00001. In all cases, as the vicinity check condition is relaxed, the tool-tilt solutions converge within fewer initial seeds.
For each test surface, the 5-axis multipoint tool position and orientation data at each tool-drop location along the all zig-zag machining passes was generated using the DSTM.
The five-axis tool positioning data computed from numerical DSTM solver was verified by plotting the entire toolpath for each of the test surfaces in a custom developed 3D toolpath simulator in MATLAB. The vector format fiveaxis tool positioning data for the toroidal tool (at each toolpath location) is defined by the location of the first point of contact S(u 1 , v 1 ), surface normal at first point of contactn p and the tool axis orientationt . Figures 6, 7, and 8 show the positioning of the toroidal tool on the test surfaces as viewed from top and bottom of the surfaces. The simulation results showed that the toroidal tool touches the surface at two points for each tool position. As discussed above, the toolpath for each surface comprises a total of 976 tool positions, but in Figs. 6b, 7b, and 8b, only every eighth tool position along each of the 16 machining passes is presented for the spin-angle value of 15 • . In the 3D plots, the bottom surface of the toroidal tool can be seen touching the three test surfaces tangentially which is visible as the circular shadow of the tool as seen in the bottom view. These circular impressions disappear as we zoom in, which confirms the gouge-free tool positioning. It is clear from the bottom views for the three test surfaces that the tool does  not protrude out of the surface in any case. Furthermore, the circular impression of the tool in the bottom views validates that the toroidal tool has a gouge-free tangential contact with the surfaces. The top views in Figs. 6a, 7a, and 8a demonstrate that the transition of tool-axis orientation along each of the machining passes is also smooth which is confirmed by the smooth orientation change of the tori along each machining pass.

Discussions
DSTM is a multipoint five-axis tool positioning method which comprised three distinct steps resulting in multiple gouge-free tool positions and orientations. Each DSTM solution has at least two gouge-free points of contact between the toroidal tool and the curved surface wherein the gap between the contact points can be altered by choosing a suitable value of the spin-angle. Three steps in the DSTM were mathematically established and solved using Maple TM [7]. The tool-drop and tool-tilt steps involve solving sets of three and four higher order rational equations respectively. In the present work, Newton's method was employed to numerically solve the simultaneous non-linear equations. A generalized method for numerically computing various partial derivative terms required in Newton's method was developed using the Richardson extrapolation technique.
The initial seeds for Newton's method were generated using a random number generator. It was observed that it is relatively easier to solve the tool-drop algorithm as it takes relatively fewer initial seeds as well as lesser number of iterations in Newton method solver as compared to the tooltilt step. This difference is attributed to the tool-drop step that require solution to only three simultaneous equations as compared to tool-tilt step that need solution to four simultaneous equations. To accelerate the computation of the second point of contact, the initial seed for the next tool position was generated by using an offset from the previous tool position. The tilt solution easily converged using such an initial seed in those cases where the surface curvature does not change significantly between two adjacent tool positions. In the event of not succeeding with such initial seeds derived from previous tool position, the random initial seeds were then used to find a solution using the numerical solver.
For a given first point of contact P , the location of the second point of contact Q in DSTM depends on the specified spin-angle. Rotating the toroidal tool by spinangle α about the surface normaln p results in a new orientation of the pseudo-insert axis given byî 3 . This results in the possibility of having a number of different gouge-free solutions for the second point of contact when the tool is tilted about the new pseudo-insert axisî 3 in the tool-tilt step. As an example, for tool-drop location T = [50, 30, 0] on the concave test surface, the first point of contact (P ) is shown as a large red color sphere in Fig. 9, while a trail of second  Fig. 9, representing the second point of contact, is influenced by the local surface curvature and effective curvature of tilted tool at first point of contact. Such a distribution is possible in relatively open surfaces. Therefore, rotation of the dropped tool through the spin-angle in the DSTM opens up the possibility of having larger side-steps of uniform cut-width between machining passes, thereby increasing the machining efficiency through reduced machining time. A gouge check algorithm [9] is used to confirm the gouge-free tangency at both points of contact.
For the tests of the DSTM algorithm for the three bicubic Bézier test surfaces, the toolpaths corresponding to different values of spin-angle (α) varying from 0 to 180 • in increments of 15 • were generated, where each toolpath comprised 976 tool positions. Figure 10 shows a subset of the tool positions (every eighth tool position) along machining passes at T x = 40 mm, 80 mm, and 120 mm for all three test surfaces. At each tool position, the first point of contact is highlighted by a sphere and the multiple solutions for the second point of contact corresponding to varying values of spin-angle are represented by an arch of polyhedra.
The shape of the distribution pattern of the second points of contact varies along the machining pass and is also dependent on the type of surface. In particular, the shape is dependent on the local surface properties. For example, in case of a concave test surface, in Fig. 10a, the arch is wide and open at the outer periphery of the surface. In this case, as the tool approaches the mid-span of the surface, the arch tends to align along the tool feed direction. This change in orientation is due to the fact that near the middle of the concave patch, the surface is relatively lower, and the tool gets the second point of contact at higher locations. Moreover, the direction of principal curvature tends to align with the feed direction as the tool proceeds towards the middle of the concave patch starting from the boundary edge. In such cases, a suitable spin-angle may help in attaining a little transverse width between the two points of contact.
Similarly, in case of the convex patch, shown in Fig. 10b, the machining pass corresponding to T x =80 mm, the trajectory of the second points of contact (at varying spinangles) is aligned along the feed direction throughout the machining pass and the trajectory of the second points of contact (around the first point of contact) has little transverse spread. This lack of transverse spread is because the direction of maximum principal curvature of the surface is aligned along the feed direction in this machining pass. This phenomenon is attributed to the fact that the convex surface patch is relatively higher in the middle than the adjacent regions. Therefore, at a given spin-angle, the tilted tool finds the gouge-free second point of contact along the feed direction only.
In case of Saddle patch, shown in Fig. 10c, the spread pattern changes its orientation as the tool progresses along a machining pass. This change in orientation is driven by changes in the local surface curvature along the machining pass.
It is observed in this study that the DSTM method is capable of achieving two gouge-free points of contact where the tool is positioned tangential to the surface. Moreover, the gap between the two tangential point of contact on the surface can be controlled by choosing a suitable value of spin-angle. Thus, using the DSTM method, a consistent cutwidth can be achieved for five-axis machining of a given parametric Bézier tensor product surface using a toroidal tool.

Conclusions
The numerical implementation of the Drop Spin and Tilt method of multipoint five-axis tool positioning for tensor product Bézier surfaces is presented in this work. The method is successfully tested on a variety of bi-quadratic, bi-cubic, bi-quartic, and bi-quintic Bézier surfaces with uniform control points. This paper highlights the unique capability of the DSTM method to generate two gouge-free points of contact between surface and toroidal tool, while the gap between two contact points can be controlled using the spin-angle.
The developed method is found to be robust as well as reliable for computing multipoint five-axis toolpath data which is demonstrated using three bi-cubic Bézier test surfaces. On the three bicubic surfaces we tested, the average number of initial seeds for our numerical method was less than 5 across all three surfaces, and the average number of Newton iterations was less than 20 except when the curvature of the tool closely matches the local surface curvature (in which case more itereations were required).
The DSTM method can be extended to other parametric surfaces with non-uniform control points. This work also opens up the avenues for future work to create efficient algorithms to develop the uniform cut width toolpaths for machining of open surfaces while simultaneously ensuring the smooth transition of tool from one tool position to the next.

Appendix A: Numerical implementation of DSTM
The DSTM method comprises three steps: tool-drop, toolspin, and tool-tilt. In the present work, Newton's method was used to determine the solutions for the tool-drop and tool-tilt steps. In the Newton's algorithm, first-and second-order partial differentials of the surface and surface normal need at each tool-drop and tool-tilt steps. To have a numerical solution that is independent of the degree of the Bèzier surface, we use the Richardson extrapolation to compute these partial differentials. The following sections discuss the implementation of the Richardson extrapolation method and Newton's method to compute the gouge free tool-drop and tool-spin solution at each tool drop location.

A.1 Richardson extrapolation for numerical differentiation
The goal of numerical differentiation is to compute an accurate estimate of the partial derivative of a function. This accuracy can be improved by decreasing the interval step size h or by using a higher order formula that employs more points. The Richardson extrapolation employs a center difference approach and offers better convergence as the truncation error drops by a factor of 4 as the step interval drops from h to h 2 . The test surfaces S used in DSTM algorithm are function of parameters u and v. Figure 11 shows the procedure developed to numerically compute the partial derivative of surface S with respect to parameter u. Similar, procedure is used to determine the partial derivative of surface S with respect to parameter v. The initial value of the interval size h step is set to 0.001 and h step is reduced by half in each subsequent iteration. The iterations continues till the progressive truncation error is less than the specified

A.2 Tool-drop algorithm
As a first step in DSTM implementation, the toroidal tool is dropped vertically through tool footprint location T (Fig. 12 gives pseudo-code for this algorithm). Equation 2 is solved by employing Newton's method to compute the location of first point of contact S(u 1 , v 1 ), and the height of the dropped tool h. Newton's method needs an initial seed to start the solution search. In the first step, a random initial solution, u, v, h], is generated and the partial derivatives required in the Newton's method are computed using the Richardson extrapolation method. Depending upon the local surface properties and the tool-drop location T , first point of contact P may lie on the outer-side or the inner-side of the cutting edge of the pseudo-insert, i.e., O 1 (see Fig. 1).
If the first point of contact P lies on the outer-side of the pseudo-insert, Eq. 3 is used to computen r . However, when the point of contact P lies on the inner side of the pseudoinsert, as shown in Fig. 13, the surface normaln p is oriented in a direction away from the tool-axist and Eq. 3 can not be used to compute the radial vectorn r . In such scenarios,n r is computed alternatively by using the following equation, Terms of the Jacobian matrix(J D ), required in the Newton's method, comprises the partial derivatives of Eq. 2 with The partial derivatives of surface S, surface normal n p , and radial vectorn r are numerically computed by employing the Richardson extrapolation method as shown in Fig. 11.
Newton's method starts with a randomly generated initial seed for {u, v, h} and attempts to converge to the solution within N I ter iterations. If Newton's method (Fig. 14) converges within the specified tolerance of Drop , then the solution is tested for gouging. If the solution yields a nongouging tool position, this is accepted as the final drop solution. If the tool is gouging with the surface or Newton's method does not converge within the specified number of iterations N I ter , then the solution is aborted and Newton' method is re-initialized with another random seed. In the gouge check step, the tool bottom is discretized into 48 (N φ ) and 36 (N θ ) divisions along φ and θ directions respectively as shown in Fig. 17.
Thus, the tool-drop algorithm yields the first nongouging point of contact P = S(u 1 , v 1 ) between the tool and the surface, the unit surface normaln p through P , the unit radial vectorn r , center of the Pseudo-Insert O 1 , and the center of the dropped toroidal tool T c .  Gouge , then the tool is considered to be gouging with surface. If the tool is found to be gouging then the solution is discarded and new trial is initiated with a new random seed. The gouge checking procedure is identical to [10] and therefore pseudo-code for the same is not given here. At any given tool-drop location T , it is also possible that a valid tool-tilt solution is not found despite considering N T ilt random seed points. In such cases, we assume that the local maximum principal curvature of the surface and the tool curvature are identical, making it infeasible to tilt the tool anymore. Therefore, in absence of a valid tilt solution, the first and second point of contact are considered to be at the same location.
Author contribution The authors contributed the work in this paper.
Funding This research was supported by NSERC.
Availability of data and material Some of the surfaces used in this work are given in the paper. The remaining surfaces and other data are available on request.
Code availability Pseudo-code for the work in this paper appears in the Appendix; the actual code is available on request.