Third-order chord error estimation for freeform contour in computer-aided manufacturing and computer numerical control systems

Abstract

The chord error employed in computer-aided manufacturing and computer numerical control systems is a crucial index to evaluate the machining accuracy of machined parts. It is usually estimated by the second-order method, that is, the osculating circle method. The second-order estimation only takes the curvature of the curve into account, which will bring about great estimation error when applying to freeform curves. In this article, a third-order method that estimates the chord error using conical helices is proposed. By investigating the geometric properties of the conical helix, it is found that there exists a conical helix that has third-order contact with the freeform curve. With the aid of this conical helix, a third-order model for estimating the chord error of freeform curves is developed. Numerical examples of three freeform curves are provided to verify the effectiveness of the proposed estimation model.

Introduction

Toolpath generation is the core technology for numerical control machining and plays an important role in guaranteeing machining accuracy and efficiency. ¹ It can be mainly classified into four categories: iso-parameter machining,^2–4 iso-plane machining,^5,6 iso-scallop machining^7–9 and constant chord error machining. ¹⁰ An essential step in all the four methods is to calculate the step length, which is generally determined by the chord error. Since the chord error describes the deviation of the actual machining surface from the desired nominal surface, it is important to calculate the chord error accurately.

Chord error is defined as the Hausdorff distance between the desired toolpath and the generated toolpath. As shown in Figure 1(a), arc $\stackrel{⌢}{\mathit{PQ}}$ is the desired toolpath and chord $\overline{\mathit{PQ}}$ is the generated toolpath. Let $\mathit{r}(s),s\in [s_i,s_{i+1}]$ denotes $\stackrel{⌢}{\mathit{PQ}}$ and $\mathit{q}(\upsilon ),\upsilon \in [0,1]$ denotes $\overline{\mathit{PQ}}$ , then chord error $\delta$ is defined as

$$\delta =max(h(\mathit{r}\left(s\right),\mathit{q}\left(\upsilon \right)),h(\mathit{q}\left(\upsilon \right),\mathit{r}\left(s\right)))$$

(1)

where $h(\mathit{r}(s),\mathit{q}(\upsilon ))=\underset{\forall s\in [s_i,s_{i+1}]}{max}\underset{\forall \upsilon \in [0,1]}{min}\Vert \mathit{r}(s)-\mathit{q}(\upsilon )\Vert$ , $h(\mathit{q}(\upsilon ),\mathit{r}(s))=\underset{\forall \upsilon \in [0,1]}{max}\underset{\forall s\in [s_i,s_{i+1}]}{min}\Vert \mathit{r}(s)-\mathit{q}(\upsilon )\Vert$ , and $\Vert ·\Vert \mathrm{stands}$ for the Euclidean norm. The chord error defined by equation (1) usually does not have an analytical expression and is difficult to calculate. Thus, rather than calculating the exact chord error, two approaches are proposed to yield approximate solutions. One is to simplify the chord error definition. Faux and Pratt ¹¹ adopted $h(\mathit{r}(s),\mathit{q}(\upsilon ))$ in equation (1) as the chord error, which has two drawbacks. First, the chord error is formulated as a nonlinear optimization problem which still does not have an analytical solution. Second, the obtained chord error is usually underestimated. ¹² The other is to locally approximate the desired contour by a simple curve and then calculate the chord error from the generated toolpath to the simple curve. The simple curve mostly used is the osculating circle.^10,13 The osculating circle has second-order contact with the desired curve $\mathit{r}(s)$ at contact point P , as depicted in Figure 1(b), where n , c and $\kappa$ denote the unit normal vector, the center of the osculating circle and the curvature of $\mathit{r}(s)$ at P , respectively. Then, the chord error is estimated as the distance from the chord $\overline{\mathit{PQ}}$ to the osculating circle, which has the following analytical expression

$$\delta =\frac{1}{\kappa }-\sqrt{\frac{1}{{\kappa }^2}-\frac{{\Vert \mathit{PQ}\Vert }^2}{4}}$$

(2)

Figure 1.

Chord error: (a) definition of chord error and (b) osculating circle estimation.

Since the osculating circle estimation has a closed-form solution, it is widely used not only in computer-aided manufacturing (CAM) system but also in computer numerical control (CNC) system. Yong and Narayanaswami ¹⁴ and Li et al. ¹⁵ applied the osculating circle estimation to check the chord error constraint at current interpolation point. Tikhon et al. ¹⁶ and Karimi and Nategh ¹⁷ adopted the osculating circle estimation to determine the maximum kinematic error permitted by the chord error tolerance. To avoid feedrate fluctuation, Lee et al. ¹⁸ and Zhao et al. ¹⁹ proposed a method to generate a jerk-limited feedrate profile. Huang and Zhu ²⁰ extended this idea to generate a jerk-continuous feedrate profile with trigonometric function. In these works, the osculating circle estimation is used twice: (1) identify critical points and (2) calculate nominal feedrates at critical points. In addition, Sun et al.^21–23 proposed several interpolation algorithms for five-axis machining, which also applied the osculating circle estimation to calculate nominal feedrates. Although the osculating circle estimation has attained extensive applications, it can result in significant estimation errors for freeform curves of which the curvatures are variable and the torsions are not zero. It can be seen from Figure 1 that for planar curves, the chord error is affected by curvature variation. For spatial curves, it is also affected by the torsion. To account for the curvature variation and torsion, a high-precision chord error estimation method using third-order approximation is proposed in this work.

The remainder of this article is organized as follows. In section “Third-order approximation of freeform curve,” the geometric properties of the conical helix are investigated, and it is proved that a freeform curve can be locally third-order approximated by a conical helix. In section “Third-order chord error estimation,” a novel third-order method for calculating chord error of freeform curve is proposed. Numerical results are given in section “Numerical examples,” and conclusions are drawn in section “Conclusion.”

Third-order approximation of freeform curve

The osculating circle method uses the osculating circle to locally approximate the freeform curve, which only takes the curvature into account. Analogously, here a conical helix is used to locally approximate the freeform curve by considering the influences of the torsion and the curvature variation.

Representation of conical helix

As shown in Figure 2(a), in coordinate frame ${\mathit{O}}_h-{\mathit{x}}_h{\mathit{y}}_h{\mathit{z}}_h$ , a conical helix ${\mathit{r}}_h(s)$ has the following parametric form

$${\mathit{r}}_h\left(s\right)={\left[\begin{array}{ccc}\beta \mathrm{as}\hfill & \alpha \mathrm{as}\hfill & s\sqrt{1-a^2}\hfill \end{array}\right]}^T$$

(3)

where auxiliary parameters $\alpha =\sin (b\ln s)/\sqrt{1+b^2}$ and $\beta =\cos (b\ln s)/\sqrt{1+b^2}$ , $a\in (0,1]$ , $b\in \Re \backslash \{0\}$ , and s stands for the arc length of the conical helix. Note that ${\mathit{r}}_h(s)$ represents a planar conical helix when $a=1$ .

Figure 2.

Two conical helices: (a) curvature variation is negative and (b) curvature variation is positive.

For any point on the conical helix, there exists a three-dimensional (3D) Frenet frame attached to it with unit tangent vector ${\mathit{t}}_h(s)$ , unit normal vector ${\mathit{n}}_h(s)$ and unit binormal vector ${\mathit{b}}_h(s)$ being its orthogonal basis. From equation (3), we have

$$\begin{array}{c}{\mathit{t}}_h\left(s\right)={\left[\begin{array}{ccc}a(\beta -\alpha b)\hfill & a(\alpha +\beta b)\hfill & \sqrt{1-a^2}\hfill \end{array}\right]}^T\hfill \\ {\mathit{n}}_h\left(s\right)={\left[\begin{array}{ccc}-\alpha -\beta b\hfill & \beta -\alpha b\hfill & 0\hfill \end{array}\right]}^T\hfill \\ {\mathit{b}}_h\left(s\right)={\left[\begin{array}{ccc}(\alpha b-\beta )\sqrt{1-a^2}\hfill & -(\alpha +\beta b)\sqrt{1-a^2}\hfill & a\hfill \end{array}\right]}^T\hfill \end{array}$$

(4)

For notation simplicity, ${\mathit{t}}_h(s)$ is denoted by ${\mathit{t}}_h$ , ${\mathit{n}}_h(s)$ by ${\mathit{n}}_h$ and ${\mathit{b}}_h(s)$ by ${\mathit{b}}_h$ in the following discussions. From equation (3), we also obtain

$${\kappa }_h\left(s\right)=\frac{\mathrm{ab}}{s},{\tau }_h\left(s\right)=\frac{b\sqrt{1-a^2}}{s},\kappa {\prime }_h\left(s\right)=-\frac{\mathrm{ab}}{s^2}$$

(5)

where ${\kappa }_h(s)$ , ${\tau }_h(s)$ and $\kappa {\prime }_h(s)$ denote the curvature, the torsion and the first-order derivative of curvature with respect to arc length s, respectively.

Since the curvature of the freeform curve is an absolute value, it can be inferred from equation (5) that parameter b is always positive. Then, the curvature variation of the conical helix defined by equation (3) is always negative, that is, $\kappa {\prime }_h(s)<0$ . If the curvature variation is positive, we can define another conical helix, which is shown in Figure 2(b), as follows

$${\mathit{r}}_h\left(s\right)={\left[\begin{array}{ccc}\beta a(S-s)\hfill & \alpha a(S-s)\hfill & (S-s)\sqrt{1-a^2}\hfill \end{array}\right]}^T$$

(6)

where $\alpha =\sin (b\ln (S-s))/\sqrt{1+b^2}$ , $\beta =\cos (b\ln (S-s))/\sqrt{1+b^2}$ and S denotes a positive constant.

Third-order approximation of freeform curve

In this section, the conical helix which has third-order contact with the freeform curve is derived.

In differential geometry, contact order is introduced to evaluate the approximation of two contact curves around the contact point. ²⁴ Two curves $\mathit{r}(s)$ and ${\mathit{r}}_h(s)$ possess contact of order k at a common point ${\mathit{r}}_0=\mathit{r}(s_0)={\mathit{r}}_h(s_h)$ if they are regular at point ${\mathit{r}}_0$ and agree there in all derivatives up to the kth order. ²⁵ The Taylor expansion of the freeform curve $\mathit{r}(s)$ at point $\mathit{r}(s_0)$ is

$$\begin{array}{l}\mathit{r}\left(s_0+\mathrm{\Delta }s\right)=\mathit{r}\left(s_0\right)+{\mathit{r}}^{\prime }\left(s_0\right)\mathrm{\Delta }s\\ +\frac{1}{2}{\mathit{r}}^{″}\left(s_0\right){\left(\mathrm{\Delta }s\right)}^2+\frac{1}{6}{\mathit{r}}^{‴}\left(s_0\right){\left(\mathrm{\Delta }s\right)}^3+o\left({\left(\mathrm{\Delta }s\right)}^3\right)\end{array}$$

(7)

The first- to the third-order derivatives of the freeform curve $\mathit{r}(s)$ with respect to arc length s at point $\mathit{r}(s_0)$ can be calculated as

$$\begin{array}{c}\mathit{r}\prime \left(s_0\right)=\mathit{t}\left(s_0\right)\hfill \\ \mathit{r}″\left(s_0\right)=\kappa \left(s_0\right)\mathit{n}\left(s_0\right)\hfill \\ \mathit{r}‴\left(s_0\right)=\kappa \prime \left(s_0\right)\mathit{n}\left(s_0\right)+\kappa \left(s_0\right)(-\kappa \left(s_0\right)\mathit{t}\left(s_0\right)+\tau \left(s_0\right)\mathit{b}\left(s_0\right))\hfill \end{array}$$

(8)

where $\mathit{t}(s_0)$ , $\mathit{n}(s_0)$ and $\mathit{b}(s_0)$ are the unit tangent vector, the unit normal vector and unit binormal vector of the freeform curve $\mathit{r}(s)$ at point $\mathit{r}(s_0)$ , respectively. For notation simplicity, they are denoted by $\mathit{t}$ , $\mathit{n}$ and $\mathit{b}$ . $\kappa (s_0)$ and $\tau (s_0)$ are the curvature and the torsion of the freeform curve $\mathit{r}(s)$ at point $\mathit{r}(s_0)$ , respectively. $\kappa \prime (s_0)$ is the first-order derivative of the curvature with respect to the arc length s at point $\mathit{r}(s_0)$ . Substituting equation (8) into equation (7) results in

$$\begin{array}{c}\multicolumn{1}{c}{\mathit{r}(s_0+\mathrm{\Delta }s)=\mathit{r}\left(s_0\right)+(1-\frac{{\kappa }^2\left(s_0\right)}{6}{\left(\mathrm{\Delta }s\right)}^2)\mathrm{\Delta }s\mathit{t}\left(s_0\right)}\\ +(\frac{\kappa \left(s_0\right)}{2}+\frac{\kappa \prime \left(s_0\right)}{6}\mathrm{\Delta }s){\left(\mathrm{\Delta }s\right)}^2\mathit{n}\left(s_0\right)\hfill \\ +\frac{\kappa \left(s_0\right)\tau \left(s_0\right)}{6}{\left(\mathrm{\Delta }s\right)}^3\mathit{b}\left(s_0\right)+o({(\mathrm{\Delta }s)}^3)\hfill \end{array}$$

(9)

According to the contact order definition and equation (9), if the torsion and the curvature variation are taken into account in locally approximating the freeform curve, the approximating curve has at least third-order contact with the freeform curve. For the third-order approximation, we have the following theorem.

Theorem 1

If the curvature $\kappa (s_0)$ of the freeform curve $\mathit{r}(s)$ at point ${\mathit{r}}_0$ satisfies $\kappa \prime (s_0)<0$ , then the freeform curve $\mathit{r}(s)$ can be locally third-order approximated by a conical helix ${\mathit{r}}_h(s)$ .

Proof

Suppose there is a conical helix ${\mathit{r}}_h\left(s\right)$ of which the Frenet frame at point ${\mathit{r}}_h\left(s_h\right)$ coincides the freeform curve $\mathit{r}\left(s\right)$ at point ${\mathit{r}}_0=\mathit{r}(s_0)$ . The transform matrix $\mathit{\Phi }$ from the Frenet frame of the conical helix ${\mathit{r}}_h(s)$ at contact point ${\mathit{r}}_0$ to that of the freeform curve $\mathit{r}(s)$ at contact point ${\mathit{r}}_0$ is

$$\mathit{\Phi }=\left[\begin{array}{ccc}\mathit{t}\hfill & \mathit{n}\hfill & \mathit{b}\hfill \end{array}\right]{\left[\begin{array}{ccc}{\mathit{t}}_h\hfill & {\mathit{n}}_h\hfill & {\mathit{b}}_h\hfill \end{array}\right]}^T$$

(10)

Then, ${\mathit{r}}_h(s_h+\mathrm{\Delta }s)-{\mathit{r}}_h(s_h)$ can be represented in the Frenet frame of freeform curve $\mathit{r}(s)$ at contact point ${\mathit{r}}_0$ as follows

$${\mathit{r}}_h\left(s_h+\mathrm{\Delta }s\right)-{\mathit{r}}_h\left(s_h\right)=\Phi {[\begin{array}{ccc}a\left({\beta }_2\left(s_h+\mathrm{\Delta }s\right)-{\beta }_1{s}_h\right)& a\left({\alpha }_2\left(s_h+\mathrm{\Delta }s\right)-{\alpha }_1{s}_h\right)& \mathrm{\Delta }s\sqrt{1-a^2}\end{array}]}^T$$

(11)

where

$$\begin{array}{c}{\alpha }_1=\sin \left(b\ln s_h\right)/\sqrt{1+b^2}\hfill \\ {\beta }_1=\cos \left(b\ln s_h\right)/\sqrt{1+b^2}\hfill \\ {\alpha }_2=\sin \left(b\ln (s_h+\mathrm{\Delta }s)\right)/\sqrt{1+b^2}\hfill \\ {\beta }_2=\cos \left(b\ln (s_h+\mathrm{\Delta }s)\right)/\sqrt{1+b^2}\hfill \end{array}$$

From equations (4) and (5), we get

$$\begin{array}{l}\left[\begin{array}{ccc}{\mathit{t}}_h& {\mathit{n}}_h& {\mathit{b}}_h\end{array}\right]=\left[\begin{array}{ccc}{{\mathit{r}}^{\prime }}_h\left(s_h\right)& {{\mathit{r}}^{″}}_h\left(s_h\right)& {{\mathit{r}}^{‴}}_h\left(s_h\right)\end{array}\right]\\ \left[\begin{array}{ccc}1& 0& \frac{a}{\sqrt{1-a^2}}\\ 0& \frac{s_h}{ab}& \frac{s_h}{a{b}^2\sqrt{1-a^2}}\\ 0& 0& \frac{s_h^2}{a{b}^2\sqrt{1-a^2}}\end{array}\right]\end{array}$$

(12)

Then, we need to prove the following equations according to equations (11) and (12)

$$\{\begin{array}{c}\mathrm{\Delta }s-\frac{a^2{b}^2}{6s_h^2}{\left(\mathrm{\Delta }s\right)}^3=\mathrm{\Delta }s-\frac{{\kappa }^2\left(s_0\right)}{6}{\left(\mathrm{\Delta }s\right)}^3+o\left({\left(\mathrm{\Delta }s\right)}^3\right)\hfill \\ \frac{\mathrm{ab}}{2s_h}{\left(\mathrm{\Delta }s\right)}^2-\frac{\mathrm{ab}}{6s_h^2}{\left(\mathrm{\Delta }s\right)}^3=\frac{\kappa \left(s_0\right)}{2}{\left(\mathrm{\Delta }s\right)}^2-\frac{\kappa \prime \left(s_0\right)}{6}{\left(\mathrm{\Delta }s\right)}^3+o\left({\left(\mathrm{\Delta }s\right)}^3\right)\hfill \\ \frac{a{b}^2\sqrt{1-a^2}}{6s_h^2}{\left(\mathrm{\Delta }s\right)}^3=\frac{\kappa \left(s_0\right)\tau \left(s_0\right)}{6}+o\left({\left(\mathrm{\Delta }s\right)}^3\right)\hfill \end{array}$$

(13)

Note that equation (13) holds if we have the following equations

$$\frac{\mathrm{ab}}{s_h}=\kappa \left(s_0\right),\frac{b\sqrt{1-a^2}}{s_h}=\tau \left(s_0\right),\frac{\mathrm{ab}}{s_h^2}=-\kappa \prime \left(s_0\right)$$

(14)

By combining equations (5) and (14), the parameters a and b of the helix conical and the helix parameter $s_h$ can be obtained

$$\begin{array}{l}a=\frac{{\kappa }_h\left(s_h\right)}{\sqrt{{\kappa }_h^2\left(s_h\right)+{\tau }_h^2\left(s_h\right)}},\\ b=-\frac{{\kappa }_h\left(s_h\right)}{{{\kappa }^{\prime }}_h\left(s_h\right)}\sqrt{{\kappa }_h^2\left(s_h\right)+{\tau }_h^2\left(s_h\right)},\\ s_h=-\frac{{\kappa }_h\left(s_h\right)}{{{\kappa }^{\prime }}_h\left(s_h\right)}\end{array}$$

(15)

Therefore, the freeform curve $\mathit{r}(s)$ can be locally third-order approximated by a conical helix determined by equation (15).

If the curvature variation of the freeform curve satisfies $\kappa \prime (s_0)<0$ , $\mathit{r}(s)$ can be locally approximated by a similar conical helix ${\mathit{r}}_h(s)$ as shown in Figure 2(b). If the curvature variation of the freeform curve satisfies $\kappa \prime (s_0)=0$ , $\mathit{r}(s)$ can be locally approximated by a cylindrical helix.

The freeform curve in Khalick’s work ²⁶ is taken as an example to illustrate the difference between the second- and third-order approximations. As can be seen in Figure 3, the conical helix that has third-order contact with the freeform curve is closer to the freeform curve than the osculating circle. Therefore, the proposed third-order approximation can attain a more accurate result in estimating the chord error.

Figure 3.

Comparison between second-order contact and third-order contact around the contact point.

Third-order chord error estimation

Since any freeform curve can be locally third-order approximated by a conical helix according to Theorem 1, the chord error of the freeform curve can be estimated by the chord error of a conical helix.

As shown in Figure 4, ${\mathit{r}}_h(s)$ is the conical helix which has third-order contact with $\mathit{r}(s)$ at ${\mathit{r}}_h(s_h)$ . The point ${\mathit{r}}_h(s_h)$ denoted by the point P is the contact point, and the point ${\mathit{r}}_h(s_h+\mathrm{\Delta }s)$ denoted by Q is the next cutter location point. $\mathit{p}(s)$ and $\mathit{q}(\upsilon )$ are, respectively, expressed as

$$\mathit{p}\left(s\right)={\mathit{r}}_h\left(s\right),s\in [s_h,s_h+\mathrm{\Delta }s]$$

(16)

$$\mathit{q}\left(\upsilon \right)=(1-\upsilon )\mathit{P}+\upsilon \mathit{Q},\upsilon \in [0,1]$$

(17)

Figure 4.

Chord error of conical helix.

According to equation (1), for a point $\mathit{p}(s)$ , when the chord error is attained, we have the following constraint for the point $\mathit{q}(\upsilon )$

$$(\mathit{p}\left(s\right)-\mathit{q}\left(\upsilon \right))·(\mathit{P}-\mathit{Q})=0$$

(18)

Solving equation (18) for parameter $\upsilon$ , we get

$$\upsilon \left(s\right)=\frac{(\mathit{p}\left(s\right)-\mathit{P})·(\mathit{Q}-\mathit{P})}{{\Vert \mathit{Q}-\mathit{P}\Vert }^2}$$

(19)

Substituting equation (19) into equation (17) yields

$$\mathit{q}\left(s\right)=\mathit{P}+(\mathit{p}\left(s\right)-\mathit{P})·\mathit{e}·\mathit{e}$$

(20)

where $\mathit{e}=(\mathit{Q}-\mathit{P})/\Vert \mathit{Q}-\mathit{P}\Vert$ .

Then, the chord error $\delta$ can be defined as

$$\delta =\underset{\forall s\in [s_h,s_h+\mathrm{\Delta }s]}{max}\sqrt{(\mathit{q}\left(s\right)-\mathit{p}\left(s\right))·(\mathit{q}\left(s\right)-\mathit{p}\left(s\right))}$$

(21)

Letting $\Psi (s)=(\mathit{q}(s)-\mathit{p}(s))·(\mathit{q}(s)-\mathit{p}(s))$ and taking its derivative with respect to arc length s, we have

$$\begin{array}{c}\multicolumn{1}{c}{\frac{d\Psi }{\mathrm{ds}}=\frac{2}{1+b^2}(1+(1-a^2)b^2)s-(1-a^2)(1+b^2)(s_h+\upsilon \mathrm{\Delta }s)}\\ \multicolumn{1}{c}{-\frac{2a^2}{1+b^2}(s_h(1-\upsilon )(\alpha {\alpha }_1+\beta {\beta }_1)+\upsilon (s_h+\mathrm{\Delta }s)(\alpha {\alpha }_2+\beta {\beta }_2))}\\ \multicolumn{1}{c}{+\frac{2a^{2}b}{1+b^2}(s_h(1-\upsilon )(\alpha {\beta }_1-{\alpha }_1\beta )+\upsilon (s_h+\mathrm{\Delta }s)(\alpha {\beta }_2-{\alpha }_2\beta ))}\end{array}$$

(22)

In the following, using Taylor expansion, we will prove that the chord error is attained when $\upsilon =1/2$ and $s=s_h+\mathrm{\Delta }s/2$ . The Taylor expansion of sine, cosine and natural logarithmic functions with respect to arc length increment $\mathrm{\Delta }s$ , up to the third-order, are used to simplify $\upsilon (s)$ and $d\Psi /\mathrm{ds}$ . $\upsilon (s)$ in equation (19) is approximately expressed as

$$\begin{array}{c}\multicolumn{1}{c}{\upsilon \left(s\right)\approx -(s-s_h)(-24(1+b^2)s_h^3{(s_h+\mathrm{\Delta }s)}^3}\\ +a^2{b}^2((-11+b^2)s^4(3s_h^2+3s_h\mathrm{\Delta }s+{\left(\mathrm{\Delta }s\right)}^2)\hfill \\ +s^2{s}_h^2((-35+b^2)s_h^2+3(-31+b^2)s_h\mathrm{\Delta }s+3(-23+b^2){\left(\mathrm{\Delta }s\right)}^2)\hfill \\ +{\mathrm{ss}}_h^3((-35+b^2)s_h^2+3(-31+b^2)s_h\mathrm{\Delta }s+3(-23+b^2){\left(\mathrm{\Delta }s\right)}^2)\hfill \\ +s^3{s}_h((79-5b^2)s_h^2-9(-15+b^2)s_h\mathrm{\Delta }s-3(-15+b^2){\left(\mathrm{\Delta }s\right)}^2)\hfill \\ +s_h^3(24s_h^3+84s_h^2\mathrm{\Delta }s+108s_h{\left(\mathrm{\Delta }s\right)}^2-(-59+b^2){\left(\mathrm{\Delta }s\right)}^3)))f^{-1}+o\left({\left(\mathrm{\Delta }s\right)}^3\right)\hfill \end{array}$$

(23)

where $\begin{array}{l}f=2s_h^3\mathrm{\Delta }s(12(1+b^2)s_h^3+36(1+b^2)s_h^2\mathrm{\Delta }s+(1+b^2)\\ (36-a^2{b}^2)s_h{(\mathrm{\Delta }s)}^2+12(1+(1-a^2)b^2){(\mathrm{\Delta }s)}^3)\end{array}$ .

And $d\Psi /\mathrm{ds}$ in equation (22) is approximately expressed as

$$\begin{array}{c}\multicolumn{1}{c}{\frac{d\Psi }{\mathrm{ds}}\approx 2s-2s_h(s_h^3+s_h^2(3+\upsilon )\mathrm{\Delta }s+3s_h(1+\upsilon ){\left(\mathrm{\Delta }s\right)}^2+(1+3\upsilon ){\left(\mathrm{\Delta }s\right)}^3){(s_h+\mathrm{\Delta }s)}^{-3}}\\ +a^2{b}^2(9s_h^4(s_h^3+s_h^2(3+\upsilon )\mathrm{\Delta }s+3s_h(1+\upsilon )\mathrm{\Delta }{s}^2+(1+3\upsilon ){\left(\mathrm{\Delta }s\right)}^3)\hfill \\ -32{\mathrm{ss}}_h^3{(s_h+\mathrm{\Delta }s)}^3+42s^2{s}_h^2{(s_h+\mathrm{\Delta }s)}^2(s_h+(1-\upsilon )\mathrm{\Delta }s)\hfill \\ -24s^3{s}_h(s_h+\mathrm{\Delta }s)(s_h^2-2s_h(-1+\upsilon )\mathrm{\Delta }s-(-1+\upsilon ){\left(\mathrm{\Delta }s\right)}^2)\hfill \\ +5s^4(s_h^3+3s_h^2(1-\upsilon )\mathrm{\Delta }s+3s_h(1-\upsilon ){\left(\mathrm{\Delta }s\right)}^2+(1-\upsilon ){\left(\mathrm{\Delta }s\right)}^3)){\left(12s_h^3{(s_h+\mathrm{\Delta }s)}^3\right)}^{-1}+o\left({\left(\mathrm{\Delta }s\right)}^3\right)\hfill \end{array}$$

(24)

Let $s=s_h+\mathrm{\Delta }s/2$ . From equation (23), we get

$${\upsilon \left(s\right)|}_{s=s_h+\mathrm{\Delta }s/2}\approx \frac{1}{2}+o\left({\left(\mathrm{\Delta }s\right)}^3\right)$$

(25)

Since the chord error in CAM and CNC systems usually falls between the millimeter and micron scale, equation (25) can be simplifies as $\upsilon =1/2$ with high accuracy. Then, from equation (24), we have

$${\frac{d\Psi }{\mathrm{ds}}|}_{s=s_h+\mathrm{\Delta }s/2,\upsilon =1/2}\approx o\left({\left(\mathrm{\Delta }s\right)}^3\right)$$

(26)

which can be simplifies as $d\Psi /\mathrm{ds}=0$ . Therefore, $s=s_h+\mathrm{\Delta }s/2$ is taken as a stationary point of $\Psi (s)$ . As shown in Figure 4, since the conical helices always lie along one side of the tangent vector, the stationary point can be regarded as the maximum point. Then, the chord error $\delta$ in equation (21) is attained when $s=s_h+\mathrm{\Delta }s/2$ and has the expression

$$\begin{array}{c}\multicolumn{1}{c}{\delta =\frac{a}{\sqrt{2(1+b^2)}}(3s_h^2+3s_h\mathrm{\Delta }s+\mathrm{\Delta }{s}^2-s_h(2s_h+\mathrm{\Delta }s)\cos \left(b\ln \left(\frac{2s_h}{2s_h+\mathrm{\Delta }s}\right)\right)}\\ {+(s_h+\mathrm{\Delta }s)(s_h\cos \left(b\ln \left(\frac{s_h}{s_h+\mathrm{\Delta }s}\right)\right)-(2s_h+\mathrm{\Delta }s)\cos \left(b\ln \left(\frac{2s_h+\mathrm{\Delta }s}{2(s_h+\mathrm{\Delta }s)}\right)\right)))}^{1/2}\hfill \end{array}$$

(27)

To analytically calculate the chord error, we should know the arc length increment $\mathrm{\Delta }s$ . Letting L denote $\overline{\mathit{PQ}}$ , we have

$$\mathrm{\Delta }{s}^2+\frac{a^2}{1+b^2}((2s_h^2+2s_h\mathrm{\Delta }s-b^2\mathrm{\Delta }{s}^2)-2s_h(s_h+\mathrm{\Delta }s)\cos \left[b\ln \frac{s_h}{s_h+\mathrm{\Delta }s}\right])-L^2=0$$

(28)

By the third-order Taylor expansion of cosine and natural logarithmic functions, we have

$$\mathrm{\Delta }{s}^2\approx L^2+o\left({\left(\mathrm{\Delta }s\right)}^3\right)$$

(29)

According to equation (29), $\mathrm{\Delta }s=L$ . Substituting equation (15) and $\mathrm{\Delta }s=L$ into equation (27), the third-order model for estimating chord error is developed.

It should be noticed that the additional computational burden using the third-order approach instead of the second-order approach mainly comes from the calculation of the third-order derivative of the freeform curve. Take a quartic non-uniform rational basis spline (NURBS) curve for example. The number of arithmetic operations introduced by the third-order derivative are 4 for addition, 13 for subtraction, 23 for multiplication and 3 for division. Hence, the increased computation will not deteriorate the real time performance.

Numerical examples

In the numerical examples, the freeform curves are first split into a series of cutter location points with a specified chord error constraint. As a result, the real chord error for each segment determined by two adjacent cutter location points is identical to the specified value. The estimated chord errors at each cutter location point are obtained through two steps: (1) approximate the freeform curve with osculating circles and conical helices and (2) estimate the chord error with a chord length equal to the corresponding segment length. As can be seen, the same cutter location point and same chord length are used for the two methods.

Example 1

In this section, a planar butterfly curve in Figure 5(a) is split into 509, 361, 292 and 253 cutter location points under the chord error constraints, that is, of 0.005, 0.01, 0.015 and 0.02 mm, respectively. The curvature of the butterfly curve is plotted in Figure 5(b). Compared with the osculating circle method, the proposed method significantly improves the estimation accuracy, as can be seen from Figures 6 and 7. The maximum improvement of estimation accuracy with the proposed method is 87.29% when the chord error constraint ${\delta }_0$ is set to be 0.005 mm, 69.92% when ${\delta }_0=0.01\mathrm{mm}$ , 90.55% when ${\delta }_0=0.015\mathrm{mm}$ and 52.23% when ${\delta }_0=0.02\mathrm{mm}$ . The reason for these significant accuracy improvements is that the osculating circle approach only takes the curvature into account and neglects the impact of the curvature variation on the chord error.

Figure 5.

Butterfly curve: (a) planar butterfly curve and (b) curvature of planar butterfly curve.

Figure 6.

Chord errors of planar butterfly curve with δ₀ = 0.005 mm: (a) real chord error and estimated chord errors provided by the osculating circle method and the proposed method, (b) estimation errors of the two methods and (c) detailed analysis of the maximum estimation error.

Figure 7.

Chord errors of planar butterfly curve—estimation errors under different chord error constraints: (a) δ₀ = 0.01 mm, (b) δ₀ = 0.015 mm and (c) δ₀ = 0.02 mm.

In addition, the estimation errors for the osculating circle method may reach to very large values even though the chord error constraints are very small. This is because the osculating circle cannot well approximate the freeform curve when the curvature variation is large. Take the cutter location point P in Figure 6 as an example, the maximum estimation error at P with the osculating circle method is 0.1919 mm, although the chord error constraint is only 0.005 mm. The osculating circle denoted by a blue dotted line and the conical helix denoted by a red solid line at P are depicted in Figure 6(c). It can be seen that the conical helix approximates the butterfly curve well, while the osculating circle bends away from the butterfly curve. As a result, the estimation error for the osculating circle method may be much larger than the specified chord error constraint. Therefore, it is essential to consider the curvature variation when estimating chord error. For the other chord error constraints, as illustrated in Figure 7(a)–(c), the maximum estimation errors with the osculating circle method are also much larger than the corresponding chord error constraints. It further demonstrates that the proposed method is more suitable for the curve with large curvature variation.

Example 2

To demonstrate the effectiveness of the proposed third-order approach for spatial curves, a flower curve in Figure 8(a) is adopted, and its curvature and torsion are plotted in Figure 8(b) and (c), respectively. Simulation results for the osculating circle method and the proposed method are plotted in Figure 9(a)–(d). Compared with the osculating circle method, the proposed method improves the estimation accuracy significantly. The maximum improvement of the estimation accuracy is 84.95% when ${\delta }_0=0.005\mathrm{mm}$ , 87.59% when ${\delta }_0=0.01\mathrm{mm}$ , 93.75% when ${\delta }_0=0.015\mathrm{mm}$ and 96.63% when ${\delta }_0=0.02\mathrm{mm}$ . Hence, the proposed method is favorable for spatial curves.

Figure 8.

Flower curve: (a) spatial flower curve, (b) curvature of spatial flower curve and (c) torsion of spatial flower curve.

Figure 9.

Chord errors of spatial flower curve—estimation errors under different chord error constraints: (a) δ₀ = 0.005 mm, (b) δ₀ = 0.01 mm, (c) δ₀ = 0.015 mm and (d) δ₀ = 0.02 mm.

In addition, the proposed method is compared with the osculating circle method in terms of computation time to further demonstrate its effectiveness for practical use. Table 1 summarizes the computation time under different chord error constraints. The computation time is the mean value of the consumed time for 1000 trials. Although the proposed method needs to calculate the third-order derivative of freeform curve, the increase in computation time is little and acceptable considering the significant improvement of estimation accuracy.

Table 1.

Comparison of computation time under different chord error constraints.

Chord error constraint (mm)	Number of cutter location points	Time (s)
		Osculating circle method	Proposed method
0.005	249	0.0291	0.0296
0.01	176	0.0206	0.0209
0.015	143	0.0169	0.0171
0.02	122	0.0144	0.0144

Example 3

A spiral toolpath used for pocket machining is adopted as an engineering example to illustrate the effectiveness of the proposed method. The spiral toolpath and its curvature are plotted in Figure 10. About 476 cutter location points are obtained with ${\delta }_0=0.005\mathrm{mm}$ . The real chord error and estimated chord errors provided by the osculating circle method and the proposed method as well as the estimation errors are plotted in Figure 11. The proposed method exhibits higher estimation accuracy than the osculating circle method. Two cutter location points where the maximum estimation errors appear are picked out for further analysis. The osculating circles and the conical helices at the two points are depicted in Figure 12. It can be seen that the conical helices approximate the spiral toolpath more accurately, which results in the better performance of the proposed method.

Figure 10.

Spiral toolpath: (a) spiral toolpath for pocket machining and (b) curvature of spiral toolpath.

Figure 11.

Chord errors of spiral toolpath with δ₀ = 0.005 mm: (a) real chord error and estimated chord errors provided by the osculating circle method and the proposed method and (b) estimation errors of the two methods.

Figure 12.

Enlarged views of cutter location points where maximum estimation errors occur.

Conclusion

Since the second-order chord error estimation only considers the influence of the curvature, it yields inaccurate estimation results for freeform curves. To consider the curvature variation and torsion, a novel third-order approach for calculating chord error of freeform curve is proposed. First, the geometric properties of the conical helix are investigated, and it is proved that a freeform curve can be third-order approximated by a conical helix. Second, with the aid of the conical helix, a third-order mathematical model for estimating the chord error of freeform curves is developed. Numerical examples for three freeform curves demonstrate that the proposed method can greatly improve the estimation accuracy in comparison with the second-order method.