This paper presents a new vector-field-based streamline smoothing method in the parametric space and a tool orientation optimization technique for five-axis machining of complex compound surfaces with torus-end cutters. Iso-planar tool path is widely used in the machining of various types of surfaces, especially for the compound surface with multiple patches, but the operations of intersecting the compound surface with a series of planes have depended considerably on the complicated optimization methods. Instead of intersecting the surface directly with planes, a novel and effective tool path smoothing method is presented, based on the iso-planar feed vector fields, for five-axis milling of a compound surface with torus-end cutters. The iso-planar feed vector field in the parametric domain is first constructed in the form of stream function that is used to generate the candidate streamlines for tool path generation. Then, a G1 blending algorithm is proposed to blend the vector fields within the adjacent parametric domains to ensure smooth transition of cross-border streamlines. Based on the smoothened streamlines in the parametric domains, pathlines along with their correspondent side sizes are selected as desirable tool paths. Concerning a high performance machining, detailed computational techniques to determine the tool axis orientation are also presented to ensure, at each cutter contact (CC) point, the torus-end cutter touches the part surface closely without gouging. Both the computational results and machined examples are demonstrated for verification and validation of the proposed methods.

## Introduction

Due to the enormous advantages in describing complex shapes, compound surfaces, which consist of multiple surface patches, have been commonly found in industrial products including consumer products, die/mold, aerospace and automotive parts, etc. Such parts are often machined on five-axis machine tools that bring some obvious benefits such as better surface finish and higher machining efficiency compared with three-axis machine tools [1,2]. When machining the surfaces with multiple patches as a whole, the iso-planar tool path method is the typical choice to eliminate the machining marks caused by cutter's entering and leaving any two adjacent patches. It is usually done by slicing the surface using a series of parallel planes. The principle of generating iso-planar tool path is simple, but the practical operations of intersecting the compound surface with a set of parallel planes are still dependent on considerably complicated optimization procedures, such as iterative intersecting line tracing [3]. On top of that, selection of correct path interval between adjacent tool paths is also challenging. Engineers typically end up selecting conservative path interval step sizes, which leads into unnecessarily long machining time [4]. Due to the challenge of computation efforts, most methods are based on three-axis machining using a simple ball-end mill cutter. However, for the situation of five-axis machining with torus-end cutters, once the iso-planar cutter contact (CC) paths are determined, sequential adjustments of tool orientations are performed to avoid local and global gouging. Such practice suffers from the lack of further consideration of resultant changes in the cutter's machining strip width (MSW) along the whole tool paths [5].

To address these problems, a novel feed vector field tool path generation method is proposed in this paper. The proposed method takes advantage of the iso-planar feed vector fields in the parametric domains for continuous and smooth tool motions and, at the same time avoids the complicated intersecting computations. Change of the MSW caused by tool orientation adjusting is also used as a feedback in the tool path generation to ensure tool path intervals between adjacent paths are optimized. In Sec. 2, literature review is provided on the related works of tool path generation and tool orientation selection for five-axis milling of compound surfaces.

## Literature Review

### Tool Path Generation for Compound Surface.

The current cutting strategies of a compound surface can be generally classified into three categories. The first feasible way is to approximate the compound surface by a set of triangular meshes, and then, the tool path generation methods, such as iso-planar methods [6,7] or the conformal mapping-based methods [811] are used to generate tool paths for milling the part surfaces. Although the tessellation of the compound surface simplifies considerably the process of tool path generation, it suffers from loss of accuracy and the degeneration of machining quality. Thus, in practical industrial applications, a compound surface can also be machined patch-by-patch [12,13]. Using this method, most of the existing tool path patterns, such as the iso-parameter tool path [14], the iso-scallop tool path [15,16], the spiral tool path [17], and the tool path along the directions of various preferred feed directions [1823], can in principle be applied to the machining of each individual patch of the given compound surface. Although it is beneficial for maintaining the machining accuracy of each patch, a difficult issue inherently in patch-by-patch machining is to find a good solution to the linkage of tool paths among different patches. Otherwise, unappealing cutting marks may be left in the neighborhood region of shared patch boundary.

To improve the machining quality in patch boundary, Cho et al. [24] and Chen and Fu [25] proposed methods to generate cross-boundary tool path with G0 position continuity based on the topological relationships among patches, so that cutter's entry/exit marks were avoided when crossing the common boundary. However, if the connection of tool paths between adjacent patches is only G0 continuity, there are possibilities of unfavorable frequent acceleration and deceleration when the cutter is moving across the patch boundary, running a risk of worse dynamic performance [26] as well as higher energy consumption [27]. To overcome the previously mentioned difficulties, in our previous work presented in Ref. [28], a vector-field-based method was proposed to generate tool path with G1 continuity between adjacent surface patches, by which the good kinematic performance and better surface quality were achieved. The result presented in our earlier work in Ref. [28] indicates the proposed method works well for compound surface with C1 continuity.

Other than the above two kinds of machining strategies on handling compound surfaces, the compound surface can be also machined as a whole with the well-known iso-planar method [29,30]. However, current iso-planar methods are mainly focused on three-axis machining. Also, as aforementioned, the current practice of finding intersecting curves between the compound surface and the sequence of intersecting planes is pretty computation intensive since it involves solving the nonlinear optimization to determine the desired CC points. Despite of this, the iso-planar tool path methods can generate continuous paths across all the adjacent patches with at least G0 continuity. When processing the compound surface consisting of multiple surface patches, the iso-planar tool path methods do have obvious advantages in generating continuous tool paths across all surface patches without interruptions of cutter lifting during machining. To maintain this advantage and simultaneously avoid the surface/plane intersecting operations, a novel streamline-type tool path generation method is proposed in this paper. The proposed vector-based tool generation method can generate continuous tool paths with G1 continuity across the boundary of two adjacent patches based on the iso-planar feed vector field as discussed in Sec. 3.

### Tool Orientation Determination.

Tool orientation is usually determined so that local and global gouging can be simultaneously avoided. In this paper, global gouging refers to the rear gouging that is beyond the neighborhood region at, or next to the CC point. On the other hand, local gouging is focused on the vicinity of the CC point, ensuring no occurrence of overcut. A common method used in local gouging avoidance is the Sturz method [31], in which a large inclination angle is set artificially to avoid the local gouging and kept constant throughout the machining process. After that, some adaptive tool orientation methods were proposed in Refs. [3234] by matching normal curvatures of the cutter's effective cutting profile with the part surface at a CC point. Further, by means of the Dupin indicatrices of the cutter's swept surface and the part surface, several new methods were developed to find the locally optimal tool orientation [35,36]. Lately, to increase the machining strip width, second-order approximation and third-order approximation methods were proposed in Refs. [37] and [38] for local tool orientation optimization. Theoretically, these local methods can produce a maximum nominal MSW without local gouging. However, since these methods ignore the possible global gouging, they might still suffer the risk of overcuts in nearby regions not immediately at, or next to the CC point.

To achieve the optimal tool orientation without global gouging, the multipoints machining method [39], rolling ball method [40], arc-intersection method [41], and penetration elimination method [42] have been developed in the last few years. These methods are focused on finding the smallest gouge-free tool inclination angle under the constraint of zero tilt angle. These results are satisfactory for global (rear) gouging avoidance but are still far from the optimal solutions in terms of the maximum MSW [43,44]. In Ref. [45], the part surface was approximated locally by a generic quadric approximation at the CC point, and then, the tool orientation angle was optimized to enlarge the MSW and to eliminate machining gouging at the same time. However, due to the inconsistency of the fitted quadric and the original surface, the results might be inaccurate in some extreme cases and still need further improvement. In Ref. [46], the tool orientations, which tend to maintain the hyper-osculating configurations when it is accessible, and smoothly switch to nearly two-contact configurations if it is not accessible due to some constraints, are generated by a global optimization procedure for five-axis machining of freeform surfaces using a flat end cylindrical tool. As a result, a good approximation quality of machining can be obtained at a relatively high computation cost. Recently, an improved rotary contact method (IRCM) [47] was proposed to adjust the tool axis by using a nested optimization of the two tool orientation angles until the cutter touches the part surface at two points. Although the method is simple and interesting, the apparent increase in computation complexity and demand of heavy computation efforts will become burdensome when trying to optimize the two orientation angles of the cutter at the same time.

In this paper, a new global tool positioning method is proposed by combining the advantages of curvature matching. Namely, the tool is initially positioned by using a local gouging avoidance method and then follows with a global gouge checking. Once the undesirable overcut occurs, the tilt angle is first optimized to decrease the penetration error. After that, adjustment of the inclination angle is performed for global gouging avoidance. This method is generally simple in implementation, and it can guarantee large MSW like other multipoint contact methods.

The remainder of this paper is arranged as follows: the vector fields fitted to the iso-planar feed directions, streamline formulation, and G1 blending of the vector fields are first presented in Sec. 3. The procedures of finding the optimal tool orientation at each CC point are discussed in detail in Sec. 4. Then, full length tool paths on the whole compound surface are generated in Sec. 5. The analytical results and experimental validations are presented in Sec. 6, followed by the concluding remarks in Sec. 7.

## Vector-Fields-Based Streamline Construction

### Vector Field Fitted to Iso-Planar Feed Directions in Parametric Domain.

This section presents the procedure of iso-planar vector field construction in parametric domains for a compound surface. As shown in Fig. 1, for a given parametric surface $S(u,v)$, there is a unique parameter direction $δu,δv$ in its parameter domain corresponding to the iso-planar feed direction at each CC point of surface. It can be analytically solved based on the differential geometry of surfaces and curves. Without the loss of generality, assume that the parametric domain of the xth patch of the compound surface has been uniformly discretized into a grid of points, and the total number of sample points is s + 1. The iso-planar feed direction $δui,δvi,i=0,1⋯s$ in the parametric plane satisfies
$Suiδui+Sviδvi=ti$
(1)
where $ti=ni×np$ is the iso-planar feed direction on the surface in the three-dimensional Cartesian space. $ni$ is the unit normal vector of the surface at the ith sample point, and $np$ is the unit normal of the intersecting planes. The iso-planar feed direction $δui,δvi$ in the parametric domain can be obtained by solving Eq. (1) as follows:
$δui=Giti·Sui−Fiti·Svi/EiGi−Fi2δvi=Eiti·Svi−Fiti·Sui/EiGi−Fi2$
(2)

where $Ei=Sui·Sui,Fi=Sui·Svi,andGi=Svi·Svi$. For the convenience of vector field construction, the calculated vector $(δui,δvi)$ by Eq. (2) needs to be normalized.

After iso-planar directions are derived in the parametric space, the vector field $φxu,v$ of the xth patch, which is described by a stream function, is constructed to capture the feed directions at an arbitrary point as
$φxu,v=∑i=0m∑j=0nNi,k(u)Nj,l(v)Ωi,jx$
(3)
where $Ni,ku$ is a k-degree B-spline basis function defined over $U$, and $Nj,lv$ is a l-degree B-spline basis function defined over $V$. $U$ and $V$ are both clamped and open knot vectors. $Ωi,jx$ is the element of the control matrix $Ωx$. Similar to the vector interpolation method presented in Ref. [48], in this paper, the streamlines are defined and derived as the level curves of a stream function. Taking the feed direction as the tangential direction of a streamline, the streamline function can be easily derived based on the fact that the gradient of a point on the stream function is orthogonal to its tangent on the level curve, as cast in the form as follows:
$δv=∂φ∂uδu=−∂φ∂v$
(4)
Thus, an overdetermined equation set for the control coefficients $Ωi,jx$ in Eq. (3) can be established by applying Eq. (4) to all the sampling points in the parameter domain, from which the unknown $Ωx$ can be solved with a least-squares method. After the stream function is determined, the correspondent streamline can be formulated as follows:
$dudv=−φvφu=ψ(u,v)$
(5)

Using Eq. (5), tool paths on the xth patch can be generated by mapping the streamlines from the parametric domain to the physical space.

### G1 Blending of Adjacent Vector Fields.

For each surface patch of a compound surface, the smoothened cross-border streamline needs to be defined. This can only be achieved when the continuity condition of adjacent vector fields is satisfied. Denote the two adjacent vector fields by $φ1(u,v)$ and $φ2(w,v)$, as shown in Fig. 2(a), which are both bi-kth degree with the shared border at u = 1 (w = 0), and the corresponding knot vectors are given by U, V and W, V. At a point on the shared border, the condition of G1 blending streamline is given as follows:
$T2=λvT1$
(6)
where T1 and T2 are tangents of streamlines for $φ1(u,v)$ and $φ2(w,v)$, respectively, and $λv$ is a ratio coefficient. For the convenience of processing, $λv$ is assumed to be a constant $λ0$. By simplifying Eq. (6) with the boundary condition, one can obtain the following equations:
$Ω0,i+12−Ω0,i2=λ0Ωm,i+11−Ωm,i1,i=0,1,…n−1$
(7a)

$Ω1,i2−Ω0,i2w¯1+k−w¯1=λ0Ωm,i1−Ωm−1,i1u¯m+k−u¯m,i=0,1,…n$
(7b)

where $w¯1,w¯1+k,u¯m$, and $u¯m+k$ are elements in the knot vectors W and U.

To make Eqs. (7a) and (7b) valid, the original control coefficients of the control matrix $Ω1$ and $Ω2$ need to be adjusted under the principle of minimum deviation. That is, the deviation between the adjusted coefficients and the original coefficients has to be minimized. In the earlier method presented in Ref. [28], the relation expressed in Eq. (7a) is processed by fixing the last row of $Ω1$ and adjusting only the first row of $Ω2$. In this paper, to handle the associated rows simultaneously, a new and more effective blending method is proposed. Let $Ω̃i,j$ be the adjusted control coefficient of $Ωi,j$. The deviation $ε1$ can be expressed as follows:
$ε1=∑i=0nΩ̃m,i1−Ωm,i12+Ω̃0,i2−Ω0,i22$
(8)
$ε1$ is minimized subject to the following condition:
$Δi2=λ0Δi1,i=0,1,…n−1$
(9)
where $Δi2$ and $Δi1$ are the unknown differences of the adjusted elements of the last row of $Ω1$ and the first row of $Ω2$, respectively. We have $Δi2=Ω̃0,i+12−Ω̃0,i2$ and $Δi1=Ω̃m,i+11−Ω̃m,i1$. Assume that $Ω̃0,02=Ω0,02$ and $Ω̃m,01=Ωm,01$, Eq. (8) can be rearranged as follows:
$ε1=∑i=1nΩm,i1−Ω̃m,01−∑j=0i−1Δj12+Ω0,i2−Ω̃0,02−λ0∑j=0i−1Δj12$
(10)
Finding the optimal solution by minimizing $ε1$ is a typical unconstrained nonlinear optimization problem, and this can be solved efficiently by some commonly used optimization techniques. The adjusted $Ω̃m,i1$ and $Ω̃0,i2$ ($i=1,2⋯n$) can be calculated as follows:
$Ω̃m,i1=Ω̃m,01+∑j=0i−1Δj1Ω̃0,i2=Ω̃0,02+λ0∑j=0i−1Δj1$
(11)
After the elements in Eq. (7a) are determined, the control coefficients involved in the m − 1th row of $Ω1$, and the second row of $Ω2$ can be further adjusted by using Eq. (7b). The adjustment deviation error $ε2$ can be defined as follows:
$ε2=Ω̃m−1,i1−Ωm−1,i12+Ω̃1,i2−Ω1,i22,i=0,1,…ns.t.:Ω̃1,i2−Ω̃0,i2w¯1+k−w¯1=λ0Ω̃m,i1−Ω̃m−1,i1u¯m+k−u¯m$
(12)
The minimum deviation error $ε2$ can be found by using the well-known Lagrange multiplier method. The final results of $Ω̃m−1,i1$ and $Ω̃1,i2$ can be explicitly formulated as follows:
$Ω̃m−1,i1=Ωm−1,i1+a2λ02Ω̃m,i1+aλ0(Ω̃0,i2−Ω1,i2)1+a2λ02Ω̃1,i2=Ω̃0,i2+a2λ02Ω1,i2+aλ0(Ω̃m,i1−Ωm−1,i1)1+a2λ02$
(13)

where $a=w¯1+k−w¯1/u¯m+k−u¯m.$ Update the original control coefficients with the adjusted ones, two vector fields can be smoothly blended so that the smoothened streamlines with G1 continuity can be achieved at the shared border of adjacent parameter domains, as shown in Fig. 2(b).

## Determination of Adaptive Tool Orientation

In this paper, the tool orientation selection method starts first with a local method by matching the normal curvature of the cutter swept surface with that of a compound surface in all directions on the tangent plane, which can be also called “comprehensive curvature matching” [35]. After that, the overcut error of the tool positioning is evaluated to see whether the global gouging appears. If the cutter does penetrate the part surface, both the two rotational axis angles are adjusted to avoid the global gouging; otherwise, the current local tool orientation is selected as the optimal one. Details of the proposed tool orientation selection method are discussed in Secs. 4.1 and 4.2.

### Tool Orientation Via Curvature Matching.

As shown in Fig. 3, a torus-end cutter moves along a tool path $c(t)$ on the surface $S(u,v)$. Denote the circular of the torus by T. The radii of T and the torus are R and r, respectively. A local coordinate system is established at the CC point, in which the tool orientation is defined by the orientation angles $(λ,ω)$. $xL and$$zL$ axes of local coordinate system are defined as the unit tangent of the tool path and the surface normal n at a CC point P; $yL$ is the cross product of $zL$ and $xL$. $λ$($0≤λ≤π/2)$ is the inclination angle measured counterclockwise around $yL$, and $ω(−π/2≤ω≤π/2)$ is the tilt angle measured counterclockwise around $zL$. To determine the geometric relationships between the cutter and the part surface, two additional coordinate systems are required, which are the tool coordinate system (TCS) and the workpiece coordinate system (WCS), respectively. By defining $xT$ as the unit vector from the origin of TCS to the CC point, $zT$ as the unit tool axis vector and $yT=zT×xT$, the swept surface $Et,θ$ of cutting torus can be mathematically expressed in the WCS by the kinematic transformation as follows:
$Et,θ=MWLMLTγ$
(14)
where
$MWL=xL1yL1zL1P1xL2yL2zL2P2xL3yL3zL3P30001,with{xL=c˙tc˙tzL=Su×SvSu×SvyL=zL×xLMLT=cos λ cos ω−sin ωsin λ cos ω−cos ωR cos λ+r sin λcos λ sin ωcos ωsin λ sin ω−sin ωR cos λ+r sin λ−sin λ0cosλR sinλ−rcos λ−10001$
$γ$ is the cutting surface depicted in the TCS, we have
$γ=R+rsin ϕcos θR+rsin ϕsin θr1−cos ϕ1T$
where $θ$($θ∈[0,2π]$) and $ϕ$($ϕ∈[0,π/2]$) describe the point location on the torus of the cutter. Since only the effective cutting curvature of the cutter swept surface at the CC point is involved, $ϕ$ is replaced with $λ$ in the following computation.
Let $καS$ denote the normal curvature of the surface in the direction $fα$ on the tangent plane, as shown in Fig. 4. $α$ is the angle from $xL$ axis to $fα$. Expression of $καS$ can be drawn from the Euler's formula shown as follows:
$καS=κmaxScos2α−β+κminSsin2α−β=11+tan2ακ0S+2τStanα+κπ/2Stan2α$
(15)

$κ0S=κmaxScos2β+κminSsin2βκπ/2S=κmaxSsin2β+κminScos2βτS=κmaxS−κminSsin βcos β$
(16)

where $κmaxS$ and $κminS$ are the principle curvatures of the part surface, $β$ specifies the maximum principle direction $fmax$ of the surface at a CC point, and $τ$ is the surface torsion in terms of $β$.

Assume that the directional curvature of the cutter swept surface $Et,θ$ is denoted by $καE$, which can be derived from the first and second partial differentiation of $Et,θ$ with respect to $θ$ and $t$. Considering Eq. (14) at the CC point, i.e., $θ=0$, one can find the following equations:
$Et=c˙tEθ=−R0xL sin ω+R0yL cos ωEtt=c¨tEtθ=R0x˙L sin ω+R0y˙L cos ω−rxL(λt cos λ sin ω+ωt sin λ cos ω)+ryLλt cos λ cos ω−ωt sin λ sin ωEθθ=−R0xL cos λ cos ω−R0yL cos λ sin ω+R0zL sin λ$
(17)
where $R0$ is $R+r sin λ$. Note that $Et$ and $Eθ$ are two vectors on the tangent plane; thus, the normal of the cutter swept surface will be identical to n as shown in the following equation:
$nE=Et×EθEt×Eθ=n$
(18)
The normal curvature of $Et,θ$ in the direction $fα$ can be calculated by
$καE=IIE/IE$
(19)
where IE and IIE are the first and second fundamental forms of $Et,θ$. Detailed formulas of finding the first and second fundamental forms, namely, IE and IIE are provided in the Appendix. Expanding Eq. (19) yields an explicit expression on $καE$
$καE=11+tan2α[κ0S+2τStanα+(sinλR+rsinλ+τSsin2ω−κ0Ssin2ω)tan2αcos2ω]$
(20)
To match the curvature of $Et,θ$ and the part surface at a CC point in an arbitrary direction $fα$, we have
$καE=καS,0≤α≤π$
(21)
The curvature matching condition can be simplified by adding Eqs. (15) and (20) into Eq. (21) as follows:
$sinλR+rsinλ+τSsin2ω−κ0Ssin2ω1cos2ω=κπ/2S$
(22)

From Eq. (22), it can be seen that a locally optimized tool orientation can be achieved by curvature matching in only one direction on the tangent plane, i.e., $α=π/2$. However, due to complexity of the surface shape, the gained tool orientation via the curvature matching method alone cannot guarantee that it is gouge-free between the cutter and the neighborhood part surface except for the vicinity of CC point. Further, global gouging avoidance procedures are required as discussed in Sec. 4.2.

### Global Gouging Detection and Avoidance.

Figure 5 depicts the geometric relation of a torus-shaped cutter and a part surface. Curve T in the WCS can be expressed as
$Tθ=MWLMLTRcosθRsinθr1T$
(23)
Imagine a set of dense enough discrete points ${Pi|i=1,2,…n}$ on the curve T are sampled, and the corresponding foot point of point $Pi$ on the part surface is represented as $Qi$, and then, the machining error at $Pi$ is defined as below:
$ei=Pi−Qi−r$
(24)

The values of $ei|i=1,2,…n$ reveal the geometrical relation clearly: a negative value implies the existence of cutter gouging, and a positive value indicates residual cusp material left on part surface after machining and zero indicates $Qi$ is a contact point. Therefore, the gouging status between the cutter and the surface can be quantified with the negative errors. The minimum negative value $Δ=min(ei|ei≤0)$ is called as the penetrating error (PE), which is a dominated factor in the adjustment of tool orientation.

Assume that the orientation angles derived from Eq. (22) is $λ0,ω0=0,$ and the cutter penetrates the part surface under the given tool orientation. The procedures of avoiding gouging are carried out as follows: First, define the minimum error for the left half ($π≤θ≤2π$) and for the right half $(0≤θ≤π)$ of T as $el$ and $er$, respectively. Generally, $el$ and $er$ are different under the initial tool orientation due to likely geometric asymmetry at the CC point as shown in Fig. 6. If $er>el$, the tilt angle should be adjusted anticlockwise, otherwise the cutter should be tilted clockwise. In the process of tilt angle rotating, the larger error will decrease while the smaller error will increase. The optimal tilt angle will lead to equal penetrating errors on both sides of the CC point, that is
$Δ=el=er$
(25)

Then, in this case, a bisection searching routine is used to find the optimal tilt angle $ωopt$ with $λ0$ fixed. To accelerate the searching speed, $el$ and $er$ are calculated quickly based on the shortest distances from the left and the right half of T to the part surface and the shortest distance between two objects can be determined efficiently with Newton iteration method. From the error distribution curves shown in Fig. 7, one can see that the smallest PE increases dramatically after optimizing the tilt angle. In addition, the approximately symmetry error distribution (Fig. 7) can be also obtained, which is beneficial to superior machining quality.

According to the definition of PE, it is clear that a smaller $Δ$ indicates smaller inclination angle adjustment to eliminate the gouging. This is also considered to be able to produce large MSW. Once the tilt angle is optimized, the PE can be taken as a function of the inclination angle $λ$ as follows:
$g=Δ(λ,ωopt,GP)$
(26)

where $GP$ denotes the point on cutter surface corresponding to the penetrating error $Δ$.

In order to completely eliminate the possible gouging, the final inclination angle is required to make sure $g≥0$. In other words, no points on the cutter will penetrate the part surface. To achieve this goal, we first retrieve the two points corresponding to $el$ and $er$ on the cutter. Denote the two points by GP1 and GP2, respectively. Then, GP1/GP2 will detach from the part surface as $λ$ increases. The optimal inclination angle $λ1$ for GP1 arises when GP1 contacts the part surface. $λ1$ can be computed with a bisection method, and the optimal inclination angle $λ2$ for GP2 can be achieved with the same process. The adjusted inclination angle $λopt$ is selected as the larger one from the calculated inclination angles ${λ1,λ2}$. Optimization of $λ$ will terminate if no gouging points are detected in the subsequent gouging check process. Otherwise, the above mentioned optimization algorithm on $λ$ will be performed iteratively until the obtained tool orientation is gouge-free.

The presented tool orientation method can guarantee at least two-point contact configuration and is similar in spirit to the multipoint contact method. Since the initial condition is properly adopted, difference between $λopt$ and $λ0$ is slight, so that the error distribution in $(λopt,ωopt)$ will be approximate to the desirable one in $λ0,ωopt,$ which denotes better surface quality after machining.

When machining the shared patch boundary, the locally optimal tool orientation may be influenced by geometric features of adjacent patches. For example, if a CC point is located at the common boundary, the principle curvatures at this point are different for the two adjacent patches. In this situation, a constant tilt angle is suggested, and the larger inclination angle computed by Eq. (22) is selected as the local tool orientation with respect to these two adjacent patches. After that, procedures of gouge checking and avoidance will be implemented to search the globally optimal tool orientation

## Tool Path Generation for Machining the Compound Surface

Streamlines in the vector field can be generated by solving Eq. (5) with a fourth-order Runge–Kutta solver. Once the initial point/seed point (SP) is determined, the Runge–Kutta integration can be implemented with a forward or backward step to trace a streamline. As long as a set of streamlines is generated with reasonably controlled intervals in the parameter domain, their corresponding pathlines on the part surface can be used as desired tool paths.

Given the machining tolerance h, two points with the position error equal to h can be found on both sides of a CC point under the optimal tool orientation. Denote their corresponding foot points on the surface with $Pl$ and $Pr$, the machining strip width L is computed by the projection distance between $Pl$ and $Pr$ along the unit vector $yL$ as follows:
$L=PlPr→·yL$
(27)

Denote surface parameters of $Pl$ and $Pr$ by $(ui,jl,vi,jl)$ and $(ui,jr,vi,jr)$, respectively, subscripts i and j denote the jth point on the ith streamline, l and r reveal their locations in terms of the feed direction at a CC point. To maximize the machining efficiency, the interval between two adjacent streamline paths should be as large as possible, meanwhile satisfying the machining tolerance. This is an iterative process as seen in Fig. 8. Assume that the streamlines are arranged in order from bottom to top in the parametric domain. The machined boundaries $Bil$, $Bir$ of the ith streamline and $Bi+1l$, $Bi+1r$ of the i + 1th streamline are formed by interpolating the points $(ui,jl,vi,jl)$, $(ui,jr,vi,jr)$, $(ui+1,jl,vi+1,jl)$, and $(ui+1,jr,vi+1,jr)$, respectively. To achieve the largest possible interval between the ith and i + 1th streamline, $Bil$ should lie on the upper side of $Bi+1r$ and be as close as possible to $Bi+1r$. This can be realized by adjusting the SP of the i + 1th streamline through trial and error. For a single patch, a set of feasible streamlines can be generated in sequence with this process until the whole parametric domain is covered.

For the streamline crossing the common boundary of adjacent patches, both the adjacent parametric domains and their associated stream functions are simultaneously taken into consideration in tool path generation. To search the SP satisfying the machining requirement of all domains crossed by the streamline, the SP for every parametric domain is first searched using the method discussed earlier; then, the intersection of a streamline with the common boundary is easily obtained based on the features of a streamline. Two intersections can be found for the boundary between two adjacent patches. The one that satisfies machining requirements of both patches is selected to be a potential SP. The final SP of the next cross-border streamline can be determined by recursively invoking the same process for every two adjacent patches until all the related domains are traversed.

Once all tool paths are generated by the above methods, the cutter location data, including the cutter center O and the tool axis vector I, can be obtained with the CC point P and the assigned orientation angles ($λ,ω)$ as expressed in the following equation:
$O=P+zLr+R/sin λ−I(r+R/tan λ)I=sin λ cos ωxL+sin λ sin ωyL+cos λzL$
(28)

## Examples

The developed tool path generation method has been coded in C++ language and implemented on a 2.8 GHz personal computer. In the first test, a compound surface with two patches as illustrated in Fig. 9(a) is used to demonstrate the validity and efficiency of the proposed method via simulation. The compound surface is modeled in Unigraphics 9.0 by stitching two Bézier patches with C1 continuity. A torus-end cutter with R = 6 mm and r = 2 mm is adopted to machine the surface under the scallop height h = 0.05 mm. Feed directions of the sample points are chosen to be parallel to the XZ plane. Then, the iso-planar feed vector fields are constructed and blended by using the method of Sec. 3. On basis of this, the streamline tool paths are generated by the proposed method and shown in Fig. 10(a). By comparison, the traditional iso-planar tool paths generated by slicing the compound surface using a series of XZ planes are shown in Fig. 10(b). It can be seen that our tool paths are similar with iso-planar tool paths in topology. The computation times listed in Table 1 show that the proposed method is more efficient than the iso-planar method, and the computation efficiency is almost doubled.

For our proposed tool path and the iso-planar tool path, the tool orientation at a CC point is determined by the proposed tool orientation method presented in Sec. 4. Then, the milling simulation of these two kinds of tool paths is performed in Unigraphics 9.0. The machined tolerance maps of these two tool paths are shown in Figs. 11(a) and 11(b), respectively. One can see that there are no overcuts on the surface. Although the machining error exceeds the required scallop height in some fractional regions, it is small enough to be neglected in the practical machining. Table 2 lists the comparison of the proposed tool orientation method with the IRCM method at ten uniformly selected points. In this test, the offset value of IRCM is set to 0, and the initial tool orientation is selected as (0, 0). Considering that it is numerically intensive when optimizing the tilt angle with a constant step as in the IRCM, a bisection method is also used for the tilt angle optimization process just as ours. Both the methods produce almost the same optimal tool orientations. However, compared with our method, the IRCM is still lower in computation efficiency due to more iterative optimizations of the two orientation angles. The experimental results, listed in Table 2, demonstrate this point. This shows that our proposed method can not only generate gouging-free tool orientation but also have higher computation efficiency.

The second compound surface, consisting of four B-spline patches as shown in Fig. 9(b), is retrieved from the computer-aided design model of a mold. On this surface, for the need of test, the traditional iso-planar tool path method is performed patch-by-patch, and the intersecting planes are selected to be parallel to the diagonal line of each patch. The generated iso-planar tool paths are shown in Fig. 12(b). Through the sampled iso-planar feed directions, the proposed method is performed, and the generated vector-field-based tool paths are shown in Fig. 12(a). Time consumption of these two kinds of tool paths is also provided in Table 1, and similar conclusion can be drawn as in Test 1. Comparing the generated tool paths, as expected, it is observed from the magnified view shown in Fig. 12(c) that the discontinuity occurs at the shared boundary for the iso-planar tool paths, while the tool paths generated by the proposed method are very smooth and continuous. An actual machining experiment is also conducted on the proposed tool path pattern and the benchmarking one. The scallop height is set to be 0.05 mm, and a torus-end cutter of radius with R = 4 mm and r = 2 mm is utilized to mill the part. The machined results are shown in Fig. 13. It reveals that the surface machined by the proposed method has better surface finish than that of the iso-planar method at the shared borders. The machining precision of the surface profile is measured by a Prismo coordinate measuring machine (CMM) in our precision machining lab. Three cross sections at X = 55 mm, Y = 45 mm, and X = 10 mm (X = 55 mm and Y = 45 mm are shared patch boundaries) are measured, respectively, on the parts shown in Figs. 13(a) and 13(b). Figure 14 shows the differences between the theoretical profile and the measured profiles. In Fig. 14, CMM data1 denote the machined profile by the proposed method, and CMM data2 are for the iso-planar method. From the results, it is seen that the maximum undercut has been constrained below the specified machining tolerance and no overcut/gouging appears. Hence, the machined parts satisfy machining requirements and validate the effectiveness of the proposed method.

## Conclusions and Future Works

In this paper, the machining of compound surface has been realized with the proposed vector-field-based smoothened streamline tool path method. This presented method takes advantages of the iso-planar feed vector and avoids the complicated intersection operations of the model surface with the intersecting planes. Since it realizes the G1 blending of feed vector fields in the parametric domains of adjacent patches, the generated tool path can cross smoothly the shared border with G1 continuity compared with the traditional iso-planar tool path. The proposed tool path method is an attempt to schedule tool paths with G1 continuity for machining the whole compound surface. The experimental results have demonstrated its validity. In addition, an effective tool positioning method is also presented by eliminating penetration errors between the cutter and the part surface. The proposed method is easy to implement, and the computed tool orientations are advantageous in achieving reasonably good machined surface quality as demonstrated in both the simulation and the actual machining experiments. Compared with IRCM method, the computation time for the optimal tool orientation is also greatly reduced.

In the future work, the emphasis should be focused on the theoretical analysis of continuous tool path blending for the arbitrary-bounded adjacent patches. We are also interested in combining the specific machine tools' kinematic constraints and the tool orientation optimization for better machining performance.

## Funding Data

• National Natural Science Foundation of China (51525501, 11290143, and 51621064).

• SCP to Dalian University of Technology (JCKY2016212A506-0102).

• National Science Foundation Grants to North Carolina State University (CMMI-1547105 and CMMI-1404916).

### Appendix: The First and Second Fundamental Formula of the Tool Swept Surface

Based on the identities in differential geometry, the first and second fundamental expressions of $Et,θ$ are given by
$IE=Et·Etdt2+2Et·Eθdtdθ+(Eθ·Eθ)dθ2IIE=Ett·ndt2+2Etθ·ndtdθ+Eθθ·ndθ2$
(A1)
where
$Et·Et=c˙t·c˙t=ISEt·Eθ=−R0sinωc˙tEθ·Eθ=R02Ett·n=c¨t·n=IISEtθ·n=−R0sinωn·c¨t/c˙t+R0cosωτEθθ·n=R0sinλ$
(A2)
IS and IIS are the first and second fundamental forms of the design surface; $dt$ and $dθ$ are parametric increments of $Et,θ$ at the CC point, which can be computed in terms of $fα$ in the tangent plane. That is
$Etdt+Eθdθ=fα$
(A3)
Given $fα=xLcos α+yLsin α$, then $dt$ and $dθ$ can be found
$dt=cosα−ω/(c˙t cos ω)dθ=sinα/R0cosω$
(A4)
By substituting Eqs. (A2)(A4) into Eq. (A1), we have
$IE=cR02IScos2ωIIE=c[R0sinλsin2αIS+R02IIScos2α−ω−2R02cos(α−ω)sinα(IISsinω−c˙tcosωτS)]$
(A5)
where $c=1/(R0c˙tcos ω)2$

## References

References
1.
Lasemi
,
A.
,
Xue
,
D.
, and
Gu
,
P.
,
2010
, “
Recent Development in CNC Machining of Freeform Surfaces: A State-of-the-Art Review
,”
Comput. Aided Des.
,
42
(
7
), pp.
641
654
.
2.
Li
,
Z.
, and
Zhu
,
L.
,
2016
, “
Mechanistic Modeling of Five-Axis Machining With a Flat End Mill Considering Bottom Edge Cutting Effect
,”
ASME J. Manuf. Sci. Eng.
,
138
(
11
), p.
111012
.
3.
Barnhill
,
R. E.
, and
Kersey
,
S. N.
,
1990
, “
A Marching Method for Parametric Surface/Surface Intersection
,”
Comput. Aided Geom. Des.
,
7
(
1–4
), pp.
257
280
.
4.
Lu
,
D.
,
Liu
,
J.
,
Zhao
,
W.
,
Lu
,
B.
,
Wu
,
D.
,
Song
,
D.
,
Xue
,
F.
, and
Cheng
,
B.
,
2017
, “
Tool Path Generation for Turbine Blades Machining With Twin Tool
,”
ASME J. Manuf. Sci. Eng.
,
139
(
11
), p.
111015
.
5.
Lee
,
Y. S.
,
1998
, “
Non-Isoparametric Tool Path Planning by Machining Strip Evaluation for 5-Axis Sculptured Surface Machining
,”
Comput. Aided Des.
,
30
(
7
), pp.
559
570
.
6.
Jun
,
C. S.
,
Kim
,
D. S.
, and
Park
,
S.
,
2002
, “
A New Curve-Based Approach to Polyhedral Machining
,”
Comput. Aided Des.
,
34
(
5
), pp.
379
389
.
7.
Hu
,
P.
,
Chen
,
L.
, and
Tang
,
K.
,
2017
, “
Efficiency-Optimal Iso-Planar Tool Path Generation for Five-Axis Finishing Machining of Freeform Surfaces
,”
Comput. Aided Des.
,
83
, pp.
33
50
.
8.
Sun
,
Y. W.
,
Guo
,
D. M.
,
Jia
,
Z. Y.
, and
Wang
,
H. X.
,
2006
, “
,”
,
28
(
7–8
), pp.
721
726
.
9.
Sun
,
Y. W.
,
Guo
,
D. M.
, and
Jia
,
Z. Y.
,
2006
, “
Spiral Cutting Operation Strategy for Machining of Sculptured Surfaces by Conformal Map Approach
,”
J. Mater. Process. Technol.
,
180
(
1–3
), pp.
74
82
.
10.
Sun
,
Y.
,
Xu
,
J.
,
Jin
,
C.
, and
Guo
,
D.
,
2016
, “
Smooth Tool Path Generation for 5-Axis Machining of Triangular Mesh Surface With Nonzero Genus
,”
Comput. Aided Des.
,
79
, pp.
60
74
.
11.
Xu
,
J.
,
Sun
,
Y.
, and
Zhang
,
L.
,
2015
, “
A Mapping-Based Approach to Eliminating Self-Intersection of Offset Paths on Mesh Surfaces for CNC Machining
,”
Comput. Aided Des.
,
62
, pp.
131
142
.
12.
Veeramani
,
D.
, and
Gau
,
Y.-S.
,
1998
, “
Models for Tool-Path Plan Optimization in Patch-by-Patch Machining
,”
Int. J. Prod. Res.
,
36
(
6
), pp.
1633
1651
.
13.
Veeramani
,
D.
, and
Gau
,
Y.
,
1998
, “
Issues in Patch-by-Patch Machining of Compound Sculptured Surfaces
,”
IIE Trans.
,
30
(
4
), pp.
341
355
.
14.
Elber
,
G.
, and
Cohen
,
E.
,
1994
, “
,”
Comput. Aided Des.
,
26
(
6
), pp.
490
496
.
15.
Suresh
,
K.
, and
Yang
,
D. C. H.
,
1994
, “
Constant Scallop-Height Machining of Free-Form Surfaces
,”
ASME J. Eng. Ind.
,
116
(
2
), pp.
253
259
.
16.
Feng
,
H. Y.
, and
Li
,
H.
,
2002
, “
Constant Scallop-Height Tool Path Generation for Three-Axis Sculptured Surface Machining
,”
Comput. Aided Des.
,
34
(
9
), pp.
647
654
.
17.
Huang
,
N.
,
Lynn
,
R.
, and
Kurfess
,
T.
,
2017
, “
Aggressive Spiral Toolpaths for Pocket Machining Based on Medial Axis Transformation
,”
ASME J. Manuf. Sci. Eng.
,
139
(
5
), p.
051011
.
18.
Kim
,
T.
, and
Sarma
,
S. E.
,
2002
, “
Toolpath Generation Along Directions of Maximum Kinematic Performance: A First Cut at Machine-Optimal Paths
,”
Comput. Aided Des.
,
34
(
6
), pp.
453
468
.
19.
Chiou
,
C. J.
, and
Lee
,
Y. S.
,
2002
, “
A Machining Potential Field Approach to Tool Path Generation for Multi-Axis Sculptured Surface Machining
,”
Comput. Aided Des.
,
34
(
5
), pp.
357
371
.
20.
Giri
,
V.
,
Bezbaruah
,
D.
,
Bubna
,
P.
, and
Choudhury
,
A. R.
,
2005
, “
Selection of Master Cutter Paths in Sculptured Surface Machining by Employing Curvature Principle
,”
Int. J. Mach. Tools Manuf.
,
45
(
10
), pp.
1202
1209
.
21.
Kumazawa
,
G. H.
,
Feng
,
H.-Y.
, and
Fard
,
M. J. B.
,
2015
, “
Preferred Feed Direction Field: A New Tool Path Generation Method for Efficient Sculptured Surface Machining
,”
Comput. Aided Des.
,
67–68
, pp.
1
12
.
22.
Ye
,
T.
,
Xiong
,
C.
,
Xiong
,
Y.
, and
Zhao
,
C.
,
2011
, “
Kinematics Constrained Five-Axis Tool Path Planning for High Material Removal Rate
,”
Sci. China Technol. Sci.
,
54
(
12
), pp.
3155
3165
.
23.
Bohez
,
E. L.
,
Makhanov
,
S. S.
, and
Munlinb
,
M.
,
2009
, “
On 5-Axis Freeform Surface Machining Optimization: Vector Field Clustering Approach
,”
,
5
(
1
), pp.
1
10
.http://www.koreascience.or.kr/article/ArticleFullRecord.jsp?cn=E1CDBZ_2005_v5n1_1&ordernum=1
24.
Cho
,
J. H.
,
Kim
,
J. W.
, and
Kim
,
K.
,
2000
, “
CNC Tool Path Planning for Multi-Patch Sculptured Surfaces
,”
Int. J. Prod. Res.
,
38
(
7
), pp.
1677
1687
.
25.
Chen
,
Z. C.
, and
Fu
,
Q.
,
2007
, “
A Practical Approach to Generating Steepest Ascent Tool-Paths for Three-Axis Finish Milling of Compound NURBS Surfaces
,”
Comput. Aided Des.
,
39
(
11
), pp.
964
974
.
26.
Sencer
,
B.
, and
Tajima
,
S.
,
2016
, “
Frequency Optimal Feed Motion Planning in Computer Numerical Controlled Machine Tools for Vibration Avoidance
,”
ASME J. Manuf. Sci. Eng.
,
139
(
1
), p.
011006
.
27.
Xu
,
K.
, and
Tang
,
K.
,
2016
, “
Optimal Workpiece Setup for Time-Efficient and Energy-Saving Five-Axis Mach-Ining of Freeform Surfaces
,”
ASME J. Manuf. Sci. Eng.
,
139
(
5
), p.
051003
.
28.
Sun
,
Y.
,
Sun
,
S.
,
Xu
,
J.
, and
Guo
,
D.
,
2017
, “
A Unified Method of Generating Tool Path Based on Multiple Vector Fields for CNC Machining of Compound NURBS Surfaces
,”
Comput. Aided Des.
,
91
, pp.
14
26
.
29.
Choi
,
B. K.
,
Lee
,
C. S.
,
Hwang
,
J. S.
, and
Jun
,
C. S.
,
1988
, “
Compound Surface Modelling and Machining
,”
Comput. Aided Des.
,
20
(
3
), pp.
127
136
.
30.
Lai
,
J. Y.
, and
Wang
,
D.-J.
,
1994
, “
A Strategy for Finish Cutting Path Generation of Compound Surfaces
,”
Comput. Ind.
,
25
(
2
), pp.
189
209
.
31.
Vickers
,
G. W.
, and
Quan
,
K. W.
,
1989
, “
Ball-Mills Versus End-Mills for Curved Surface Machining
,”
ASME J. Eng. Ind.
,
111
(
1
), pp.
22
26
.
32.
Jensen
,
C. G.
, and
Anderson
,
D. C.
,
1993
, “
Accurate Tool Placement and Orientation for Finished Surface Machining
,”
J. Des. Manuf.
,
3
(
4
), pp.
251
261
.
33.
Jensen
,
C. G.
,
Anderson
,
D. C.
, and
Mullins
,
S. H.
,
1993
, “
Scallop Elimination Based on Precise 5-Axis Tool Placement, Orientation, and Step-Over Calculations
,”
,
65
(
2
), pp.
535
544
.
34.
Rao
,
N.
,
Bedi
,
S.
, and
Buchal
,
R.
,
1996
, “
Implementation of the Principal-Axis Method for Machining of Complex Surfaces
,”
,
11
(
4
), pp.
249
257
.
35.
Rao
,
A.
, and
Sarma
,
R.
,
2000
, “
On Local Gouging in Five-Axis Sculptured Surface Machining Using Flat-End Tools
,”
Comput. Aided Des.
,
32
(
7
), pp.
409
420
.
36.
Yoon
,
J.-H.
,
Pottmann
,
H.
, and
Lee
,
Y.-S.
,
2003
, “
Locally Optimal Cutting Positions for 5-Axis Sculptured Surface Machining
,”
Comput. Aided Des.
,
35
(
1
), pp.
69
81
.
37.
Gong
,
H.
,
Fang
,
F. Z.
,
Hu
,
X. T.
,
Cao
,
L.-X.
, and
Liu
,
J.
,
2010
, “
Optimization of Tool Positions Locally Based on the BCELTP for 5-Axis Machining of Free-Form Surfaces
,”
Comput. Aided Des.
,
42
(
6
), pp.
558
570
.
38.
Zhu
,
L.
,
Ding
,
H.
, and
Xiong
,
Y.
,
2010
, “
Third-Order Point Contact Approach for Five-Axis Sculptured Surface Machining Using Non-Ball-End Tools—Part II: Tool Positioning Strategy
,”
Sci. China Technol. Sci.
,
53
(
8
), pp.
2190
2197
.
39.
Warkentin
,
A.
,
Ismail
,
F.
, and
Bedi
,
S.
,
1998
, “
Intersection Approach to Multi-Point Machining of Sculptured Surfaces
,”
Comput. Aided Geom. Des.
,
15
(
6
), pp.
567
584
.
40.
Gray
,
P. J.
,
Bedi
,
S.
, and
Ismail
,
F.
,
2003
, “
Rolling Ball Method for 5-Axis Sur-Face Machining
,”
Comput. Aided Des.
,
35
(
4
), pp.
347
357
.
41.
Gray
,
P. J.
,
Bedi
,
S.
, and
Ismail
,
F.
,
2005
, “
Arc-Intersect Method for 5-Axis Tool Positioning
,”
Comput. Aided Des.
,
37
(
7
), pp.
663
674
.
42.
Hosseinkhani
,
Y.
,
Akbari
,
J.
, and
Vafaeesefat
,
A.
,
2007
, “
Penetration–Elimination Method for Five-Axis CNC Machining of Sculptured Surfaces
,”
Int. J. Mach. Tools Manuf.
,
47
(
10
), pp.
1625
1635
.
43.
Barakchi Fard
,
M. J.
, and
Feng
,
H.-Y.
,
2009
, “
Effect of Tool Tilt Angle on Machining Strip Width in Five-Axis Flat-End Milling of Free-Form Surfaces
,”
,
44
(
3–4
), pp.
211
222
.
44.
Barakchi Fard
,
M. J.
, and
Feng
,
H.-Y.
,
2010
, “
Effective Determination of Feed Direction and Tool Orientation in Five-Axis Flat-End Milling
,”
ASME J. Manuf. Sci. Eng.
,
132
(
6
), p.
061010
.
45.
Fan
,
J.
, and
Ball
,
A.
,
2014
, “
Flat-End Cutter Orientation on a Quadric in Five-Axis Machining
,”
Comput. Aided Des.
,
53
, pp.
126
138
.
46.
Kim
,
Y.-J.
,
Elber
,
G.
,
Bartoň
,
M.
, and
Pottmann
,
H.
,
2015
, “
Precise Gouging-Free Tool Orientations for 5-Axis CNC Machining
,”
Comput. Aided Des.
,
35
, pp.
220
229
.
47.
Fan
,
W.
,
Ye
,
P.
,
Fang
,
C.
, and
Shi
,
S.
,
2013
, “
Improved Rotary Contact Method for 5-Axis Sculptured Surfaces Machining
,”
ASME J. Manuf. Sci. Eng.
,
135
(
5
), p.
051002
.
48.
Preusser
,
A.
,
1986
, “
Computing Area Filling Contours for Surfaces Defined by Piecewise Polynomials
,”
Comput. Aided Geom. Des.
,
3
(
4
), pp.
267
279
.