Jerk Limited Feedrate Profile Optimization for Nurbs Toolpaths in Cnc Machines with H-Bot Kinematics

ERWIŃSKI, Krystian; USTASZEWSKI, Patryk

doi:10.2478/ama-2025-0075

Full Article

1.

INTRODUCTION

Computerized Numerically Controlled CNC machines are the backbone of the manufacturing industry. There are many types of such machines, including milling machines, lathes, laser cutters and industrial 3D printers. Motion of the end effector (milling spindle, laser, print head), which affects the workpiece, is the result of the coordinated motion of individual machine axes. Each axis is driven by one or more electric motors controlled directly by the computerized numerical controller.

Gantry type CNC machines are usually used as milling plotters or cutters using laser, plasma or oscillating tangential knife. The end effector, such as a milling spindle is mounted on a horizontal beam mounted on two supports which move along the length of the workpiece table, called a portal. The beam also mounts the transverse axis which allows for lateral motion of the end effector. These types of machines are usually used for cutting or milling large sheets of metal, plastic, leather. Compared to traditional milling machines the gantry type offers large work space and higher feedrates.

Typically the gantry CNC is constructed as a machine with cartesian kinematics. Due to large work area, the beam is driven by two parallel drivetrains, usually based on a rack and pinion mechanism. The lateral axis also features a rack and pinion or ballscrew drive train depending on its length. The downside of this approach is that the axes have to drive both the mass of the beam, the rack and pinion mechanism as well as the motors with their gearboxes. The increased mass conflicts with the requirement of achieving high speeds and accelerations.

An alternative solution for a gantry CNC machine is the H-Bot gantry kinematic configuration. It has functionally the same applications – mainly cutting of sheet materials using laser, plasma, tangential knife or 2D milling. Both types of machines feature a portal beam riding on two parallel guideways on which a transverse axis is mounted. In the H-Bot machine the axes are actuated by a single toothed belt guided by a system of pulleys in the form of a letter H (hence the name). The belt is driven by two servo motors usually with reduction gearboxes. Contrary to the traditional cartesian gantries these motors are mounted on the machine frame and do not move with the axes. The axes do not require mounting of heavy rack and pinion mechanism and motor cable guides. This significantly decreases the weight driven by the machine motors. This can improve the machine dynamics as higher accelerations and velocities can be achieved with the same power of the motors. Compared to rack and pinion mechanism the toothed belt drive train features much less backlash when properly tensioned. Furthermore the H-Bot gantry requires only two motors, servo drives and gearboxes instead of three which are used in cartesian gantry CNC machines. In practice the H-bot gantry can replace the cartesian gantry especially when used in high speed cutting applications due to its lighter weight, higher dynamic capabilities and lower costs. A schematic comparing the cartesian and H-bot gantries is shown in (Fig. 1).

Feedrate profile defines the feedrate (end-effector velocity tangent to the toolpath) as a function of time or the toolpath parameter. Feedrate profile generation or feedrate planning is performed by taking into account limitations of the machines capabilities such as drive power, maximum drive speed, drivetrain parameters and machine dynamics [1].

In its simplest form, the feedrate profile limits only tangent velocity and acceleration, and the profile shape is a trapezoid [2]. In high-performance machines, the feedrate profile also limits jerk, which is the time derivative of acceleration. The feedrate profile is then represented by an S-Curve profile [3]. Jerk limitation causes the acceleration to change gradually, which in turn limits the rate of change of drive torque.

Many commercially available CNC machines use linear or circular segments for toolpath definition. Complex shapes are usually defined as a large number of short linear segments. Due to discontinuities, motion along such paths would result in discontinuous feedrate and acceleration, which in turn could cause vibration and large following errors. To avoid this, the paths are usually smoothed internally by the controller. Alternatively, in state-of-the-art solutions, polynomial paths are used as toolpath definitions. These are usually Non-Uniform Rational B-Splines (NURBS) as they are also used in CAD/CAM software to describe complex shapes. Toolpath definition as polynomial NURBS curves allows dynamic and smooth motion along complex paths.

To achieve high performance, feedrate profile generation with just feedrate, tangent acceleration and jerk limits may not be sufficient. When using NURBS toolpaths, due to their non-linear nature, the acceleration and jerk of the machine’s axes may be violated even if tangent limits are kept. The problem of maximizing the feedrate while avoiding axial velocity, acceleration and jerk limit violations for all axes is a difficult problem, especially if the toolpath has highly variable curvature. This difficulty is increased if the machine kinematics is non-Cartesian due to an additional relationship between Cartesian and drive axes.

The solutions can be broadly divided into two categories. The indirect approach usually defines a Feedrate Limit Curve or Feedrate Limit Function, which combines axial velocity limits and linearised axial acceleration and jerk limits into a single feedrate limit. Then the feedrate profile is planned under the Feedrate Limit Curve, usually using S-Curve segments [4][9][19]. This approach is usually computationally efficient but may not achieve maximum performance or may violate axial limits. Some researchers suggest using a higher-order jerk continuous feedrate profile, which should further improve motion smoothness, especially for machines with flexible drive trains [20][21].

The direct approach solves the optimization problem using a constrained optimization algorithm. The feedrate profile is adjusted by the optimization algorithm to maximize the feedrate while not violating the axial velocity, acceleration and jerk limits [10][11][22]. This approach can usually achieve higher feedrates, which implies shorter manufacturing times, but at the cost of higher computational costs. The optimization approach is also much more capable in constraining axial velocity, acceleration and jerk.

In many cases the proposed methods consider jerk-limited feedrate profiles in which the jerk may be discontinuous. For light-weight low-stiffness machines and those which feature elasticity in the drive train, like the belt-driven H-bot, a higher order jerk-continuous feedrate profile may be preferable. In [20] a trigonometric-function-based jerk-continuous feedrate optimization method is proposed. In [23] a quartic feedrate profile is used. In [24] a snap limited linear feedrate optimization is used for an articulated robot. These methods do not consider axis limits or include them indirectly. Direct inclusion of axial limits and optimization of a jerk-continuous feedrate profile is usually a complex optimization problem usually requiring significant computational effort.

In most cases, the feedrate planning problem was investigated for 3-axis or 5-axis Cartesian machines. In this paper, the optimization approach is used to generate a jerk-continuous cubic feedrate profile considering axial limitations for a non-Cartesian H-Bot CNC machine with toolpaths defined as cubic NURBS curves. To improve feedrate generation computation, the feedrate profile is initially generated as an S-Curve and then converted to a cubic B-Spline representation. Optimization of the B-Spline is then performed so that the axial velocity, acceleration and jerk of the H-Bot machine are continuous and their limits are not violated. The Particle Swarm Optimization algorithm is used to optimize the shape of the profile, and the Augmented Lagrangian Method is used to handle axial constraints.

2.

NURBS TOOLPATH INTERPOLATION FOR H-BOT

Non-Uniform Rational B-Splines (NURBS) are polynomial splines defined as: 1 $C (u) = \frac{\sum_{i = 1}^{k} N_{i,n} (u) ω_{i} P_{i}}{\sum_{i = 1}^{k} N_{i,n} (u) ω_{i}}$ {\bf{C}}({\bf{u}}) = {{\sum\nolimits_{{\rm{i}} = 1}^{\rm{k}} {{{\rm{N}}_{{\rm{i,n}}}}({\rm{u}}){{\rm{\omega }}_{\rm{i}}}\;{{\bf{P}}_{\rm{i}}}} } \over {\sum\nolimits_{{\rm{i}} = 1}^{\rm{k}} {{{\rm{N}}_{{\rm{i,n}}}}({\rm{u}}){{\rm{\omega }}_{\rm{i}}}} }} where: C(u) – point on the NURBS curve, N_i,n(u)– basis functions, u – curve parameter, P_i – control point, ω_i – control point weight, k – number of control points

The control points and their weights define the shape of the curve. An example of a NURBS toolpath with marked control points is shown in (Fig. 2). The unitless parameter u defines a unique position on the curve. By increasing the parameter from 0 (start of curve) to 1 (end of curve), the curve is traversed. For defining a toolpath of a CNC machine, a cubic curve is usually used, which has C2 parametric continuity, which means that the first and second derivatives of the curve are continuous. This guarantees that during motion along the curve, the feedrate and acceleration are continuous, and the tangent jerk is finite. In practice, evaluation of the NURBS curve is performed using the numerically stable de Boor’s algorithm [6]. This method is also used to obtain the derivatives of the curve with respect to the curve parameter.

Interpolation of the NURBS toolpath in CNC control systems is the process of determining consecutive positions C(u_i) on the curve by determining consecutive parameter values in constant time intervals. The increment of the curve parameter Δu_i should correspond to an increment of the curve arc-length Δs_i which in turn corresponds to the current value of feedrate V_t(u_i). This feedrate is determined by the feedrate profile. An illustration of this process is presented in (Fig. 3).

The relationship between the parameter increment and arc length is non-linear and unique for each curve. During each interpolation step, this relationship has to be approximated. There are two main approaches to this problem. The first approach is to use a Taylor series expansion to determine consecutive values of the parameter u. This has the advantage of using a single closed-form formula and low computational complexity. The Taylor series method has the disadvantage of introducing feedrate fluctuations and interpolation errors if the feedrate is large and the curvature of the NURBS toolpath has high variations, such as sharp corners [8].

The second approach is a predictor-corrector method in which the initial parameter increment is verified by performing a numerical integration of the NURBS curve. The increment is then corrected so that the difference between the initial arc-length increment and the desired increment is minimised. This allows for interpolation of complex shapes with higher accuracy and better conformity to the desired feedrate profile [5]. In this paper, a second-order Runge-Kutta predictor-corrector (RK2-PC) method is used [7].

Interpolated Cartesian points on the curve are converted to axial position commands of the H-Bot using the inverse kinematics transformation: 2 $[\begin{array}{l} q_{1} \\ q_{2} \end{array}] = \frac{G}{R} [\begin{matrix} - 1 & 1 \\ 1 & 1 \end{matrix}] \cdot [\begin{array}{l} x \\ y \end{array}]$ \left[ {\matrix{ {{q_1}} \hfill \cr {{q_2}} \hfill \cr } } \right] = {G \over R}\left[ {\matrix{ { - 1} & 1 \cr 1 & 1 \cr } } \right]\;\cdot\;\left[ {\matrix{ x \hfill \cr y \hfill \cr } } \right] where: R – radius of the drive pulley in mm, G – motor gear ratio, q₁, q₂ – position of the drive axes, x,y – Cartesian position of the end effector.

In order to perform feedrate planning with axial limts, axial velocity, acceleration and jerk have to be determined based on feedrate, tangential acceleration and tangential jerk: 3 $V_{c} = [\begin{matrix} \dot{x} \\ \dot{y} \end{matrix}] = [\begin{array}{l} v_{x} \\ v_{y} \end{array}] = C_{s} \cdot V_{t}$ {{\bf{V}}_{\bf{c}}} = \left[ {\matrix{ {\dot x} \cr {\dot y} \cr } } \right] = \left[ {\matrix{ {{v_x}} \hfill \cr {{v_y}} \hfill \cr } } \right] = {{\bf{C}}_{\bf{s}}}\;\cdot\;{V_t} 4 $A_{c} = [\begin{array}{l} a_{x} \\ a_{y} \end{array}] = C_{ss} \cdot V_{t}^{2} + C_{s} \cdot A_{t}$ {{\bf{A}}_{\bf{c}}} = \left[ {\matrix{ {{a_x}} \hfill \cr {{a_y}} \hfill \cr } } \right] = {{\rm{C}}_{{\rm{ss}}}}\;\cdot\;V_t^2 + {{\rm{C}}_{\rm{s}}}\;\cdot\;{A_t} 5 $J_{c} = [\begin{array}{l} j_{x} \\ j_{y} \end{array}] = C_{sss} \cdot V_{t}^{3} + 3 C_{ss} \cdot V_{t} A_{t} + C_{s} J_{t}$ {{\bf{J}}_{\bf{c}}} = \left[ {\matrix{ {{j_x}} \hfill \cr {{j_y}} \hfill \cr } } \right] = {{\rm{C}}_{{\rm{sss}}}}\;\cdot\;V_t^3 + 3{{\rm{C}}_{{\rm{ss}}}}\;\cdot\;{V_t}{A_t} + {{\rm{C}}_{\rm{s}}}{{\rm{J}}_{\rm{t}}} where: C_s, C_ss, C_sss are first, second and third derivative of the NURBS curve with respect to the curve arc length,

V_c, A_c, J_c are vectors of velocity, acceleration and jerk along Cartesian axes. To determine toolpath derivatives with respect to arclength the following formulas are used: 6 $C_{s} = C^{'} u_{s}$ {{\bf{C}}_{\bf{s}}} = {{\bf{C}}^\prime }{u_s} 7 $C_{s s} = C^{''} u_{s}^{2} + C^{'} u_{ss}$ {{\bf{C}}_{{\bf{ss}}}} = {{\rm{C}}^{\prime \prime }}{\rm{u}}_s^2 + {{\rm{C}}^\prime }{{\rm{u}}_{{\rm{ss}}}} 8 $C_{s s s} = C^{'''} u_{s}^{3} + 3 C^{''} u_{ss} u_{s} + C^{'} u_{sss}$ {{\bf{C}}_{{\bf{sss}}}} = {{\bf{C}}^{\prime \prime \prime }}{\rm{u}}_s^3 + 3{{\bf{C}}^{\prime \prime }}{{\rm{u}}_{{\rm{ss}}}}{{\rm{u}}_{\rm{s}}} + {{\bf{C}}^\prime }{{\rm{u}}_{{\rm{sss}}}} where: u_s, u_ss, u_sss – first, second and third derivatives of the toolpath curve parameter with respect to the arc-length.

C′, C″, C″′ – first, second and third derivatives of the NURBS curve with respect to the curve parameter.

Parameter derivatives with respect to arc-length are defined as: 9 $u_{s} = \frac{1}{| C^{'} |}$ {u_s} = {1 \over {\left| {{C^\prime }} \right|}} 10 $u_{s s} = - C^{'} \cdot C^{''} u_{s}^{4}$ {u_{ss}} = - {C^\prime }\;\cdot\;{C^{\prime \prime }}u_s^4 11 $u_{s s s} = 4 {(C^{'} \cdot C^{''})}^{2} u_{s}^{7} - (C^{'' 2} + C^{'} \cdot C^{'''}) u_{s}^{5}$ {u_{sss}} = 4{\left({{C^{\prime \;}}\cdot\;{C^{\prime \prime }}} \right)^2}u_s^7 - \left({{C^{\prime \prime 2}} + {C^\prime }\cdot\;{C^{\prime \prime \prime }}} \right)u_s^5

Axial velocity, acceleration and jerk can be expressed as: 12 $V_{q} = [\begin{array}{l} v_{1} \\ v_{2} \end{array}] = 𝓙^{- 1} \cdot V_{c}$ {{\bf{V}}_{\bf{q}}} = \left[ {\matrix{ {{v_1}} \hfill \cr {{v_2}} \hfill \cr } } \right] = {{\bf{J}}^{{\bf{ - 1}}}}\cdot\;{{\bf{V}}_{\bf{c}}} 13 $A_{q} = [\begin{array}{l} a_{1} \\ a_{2} \end{array}] = 𝓙^{- 1} \cdot A_{c} + \dot{𝓙^{- 1}} \cdot V_{c}$ {{\bf{A}}_{\bf{q}}} = \left[ {\matrix{ {{a_1}} \hfill \cr {{a_2}} \hfill \cr } } \right] = {{\bf{J}}^{ - {\bf{1}}}}\cdot\;{{\bf{A}}_{\bf{c}}} + \mathop {{{\bf{J}}^{ - {\bf{1}}}}}\limits^{\bf{.}} \cdot\;{{\bf{V}}_{\bf{c}}} 14 $J_{q} = [\begin{array}{l} j_{1} \\ j_{2} \end{array}] = 𝓙^{- 1} \cdot J_{c} + 2 \dot{𝓙^{- 1}} \cdot A_{c} + \overset{..}{𝓙^{- 1}} \cdot V_{c}$ {{\bf{J}}_{\bf{q}}} = \left[ {\matrix{ {{j_1}} \hfill \cr {{j_2}} \hfill \cr } } \right] = {{\bf{J}}^{{\bf{ - 1}}}}\cdot{{\bf{J}}_{\bf{c}}} + 2\mathop {{{\bf{J}}^{ - {\bf{1}}}}}\limits^{\bf{.}} \cdot{{\bf{A}}_{\bf{c}}} + \mathop {{{\bf{J}}^{ - {\bf{1}}}}}\limits^{..} \cdot\;{{\bf{V}}_{\bf{c}}} where: 𝓙⁻¹ is the inverse Jacobian matrix and $\dot{𝓙^{- 1}}, \overset{..}{𝓙^{- 1}}$ \mathop {{{\bf{J}}^{ - {\bf{1}}}}}\limits^{\bf{.}} ,\;\mathop {{{\bf{J}}^{ - {\bf{1}}}}}\limits^{..} are its first and second time derivatives. The Jacobian is the matrix of partial derivatives of axial positions q₁, q₂ with respect to Cartesian positions x,y.

The final formulas that of axial velocities, accelerations and jerks are expressed as: 15 $V_{q} = [\begin{array}{l} v_{1} \\ v_{2} \end{array}] = \frac{G}{R} [\begin{matrix} - 1 & 1 \\ 1 & 1 \end{matrix}] \cdot [\begin{array}{l} v_{x} \\ v_{y} \end{array}]$ {{\bf{V}}_{\bf{q}}} = \left[ {\matrix{ {{v_1}} \hfill \cr {{v_2}} \hfill \cr } } \right] = {G \over R}\left[ {\matrix{ { - 1} & 1 \cr 1 & 1 \cr } } \right]\;\cdot\;\left[ {\matrix{ {{v_x}} \hfill \cr {{v_y}} \hfill \cr } } \right] 16 $A_{q} = [\begin{array}{l} a_{1} \\ a_{2} \end{array}] = \frac{G}{R} [\begin{matrix} - 1 & 1 \\ 1 & 1 \end{matrix}] \cdot [\begin{array}{l} a_{x} \\ a_{y} \end{array}]$ {{\bf{A}}_{\bf{q}}} = \left[ {\matrix{ {{a_1}} \hfill \cr {{a_2}} \hfill \cr } } \right] = {G \over R}\left[ {\matrix{ { - 1} & 1 \cr 1 & 1 \cr } } \right]\;\cdot\;\left[ {\matrix{ {{a_x}} \hfill \cr {{a_y}} \hfill \cr } } \right] 17 $J_{q} = [\begin{array}{l} j_{1} \\ j_{2} \end{array}] = \frac{G}{R} [\begin{matrix} - 1 & 1 \\ 1 & 1 \end{matrix}] \cdot [\begin{array}{l} j_{x} \\ j_{y} \end{array}]$ {{\bf{J}}_{\bf{q}}} = \left[ {\matrix{ {{j_1}} \hfill \cr {{j_2}} \hfill \cr } } \right] = {G \over R}\left[ {\matrix{ { - 1} & 1 \cr 1 & 1 \cr } } \right]\;\cdot\;\left[ {\matrix{ {{j_x}} \hfill \cr {{j_y}} \hfill \cr } } \right]

The formulas described above are used to compute axial kinematic limit violations during feedrate profile optimization.

3.

FEEDRATE OPTIMIZATION ALGORITHM

Feedrate optimization is performed in two phases. In the first phase, a jerk-limited S-Curve feedrate profile is generated with maximum allowable feedrate, tangential acceleration and tangential jerk. The initial polynomial profile is then converted to a B-Spline form defined by a set of control points, similar to the NURBS toolpath. Optimization is performed, using the Particle Swarm Optimization (PSO) algorithm, which maximizes the feedrate while taking into account constraints on axial velocity, acceleration and jerk. A block schematic of the feedrate optimization process is presented in (Fig. 4).

3.1.

S-Curve initial guess

S-curve feedrate profile, presented in (Fig. 5) consists of 7 second-order polynomial segments with trapezoidal acceleration profile and discontinuous jerk profile.

The profile can be mathematically defined as a function of time: 18 $V (t) = {\begin{array}{l} J_{m a x} t^{2} + V_{s t a r t} & 0 \geq t \geq τ_{1} \\ A_{1} t + V_{1} & τ_{1} > t \geq τ_{2} \\ - \frac{1}{2} J_{m a x} t^{2} + A_{2} t + V_{2} & τ_{2} > t \geq τ_{3} \\ V_{m a x} & τ_{3} \geq t > τ_{4} \\ - \frac{1}{2} J_{m a x} t^{2} + V_{4} & τ_{4} > t \geq τ_{5} \\ A_{5} t + V_{5} & τ_{5} > t \geq τ_{6} \\ \frac{1}{2} J_{m a x} t^{2} + A_{6} t + V_{6} & τ_{6} > t \geq τ_{7} \end{array}$ V(t) = \left\{ {\matrix{ {{J_{max}}{t^2} + {V_{start }}} \hfill & {0 \ge t \ge {\tau _1}} \hfill \cr {{A_1}t + {V_1}} \hfill & {{\tau _1} > t \ge {\tau _2}} \hfill \cr { - {1 \over 2}{J_{max}}{t^2} + {A_2}t + {V_2}} \hfill & {{\tau _2} > t \ge {\tau _3}} \hfill \cr {{V_{max}}} \hfill & {{\tau _3} \ge t > {\tau _4}} \hfill \cr { - {1 \over 2}{J_{max}}{t^2} + {V_4}} \hfill & {{\tau _4} > t \ge {\tau _5}} \hfill \cr {{A_5}t + {V_5}} \hfill & {{\tau _5} > t \ge {\tau _6}} \hfill \cr {{1 \over 2}{J_{max}}{t^2} + {A_6}t + {V_6}} \hfill & {{\tau _6} > t \ge {\tau _7}} \hfill \cr } } \right. where: V_max, A_max, J_max – maximum feedrate, tangential acceleration and jerk, V_1-6 – feedrate at the end of each phase, V_start – initial feedrate, A_1-6 – acceleration at the end of each phase, T_1–7 – time at the end of each phase

To plan an S-curve feedrate profile on a certain distance Sseg, the distance traversed during each phase is first computed assuming maximum feedrate, acceleration and jerk. The distance in the constant feedrate phase is adjusted so that the sum of distances of each phase matches the desired displacement.19 ${\begin{array}{l} S_{1} = \frac{1}{6} J_{m a x} T_{1}^{3} + V_{s t a r t} T_{1} \\ S_{2} = \frac{1}{2} A_{1} T_{2}^{2} + F_{1} T_{2} \\ S_{3} = - \frac{1}{6} J_{m a x} T_{3}^{3} + \frac{1}{2} A_{2} T_{3}^{2} + F_{2} T_{3} \\ S_{5} = - \frac{1}{6} J_{m a x} T_{5}^{3} + V_{m a x} T_{5} \\ S_{6} = - \frac{1}{2} A_{5} T_{6}^{2} + V_{5} T_{6} \\ S_{7} = \frac{1}{6} J_{m a x} T_{7}^{3} - \frac{1}{2} A_{6} T_{7}^{2} + V_{6} T_{7} \\ S_{4} = S_{seg} - (S_{1} + S_{2} + S_{3} + S_{5} + S_{6} + S_{7}) \end{array}$ \left\{ {\matrix{ {{S_1} = {1 \over 6}{J_{max}}T_1^3 + {V_{start{\rm{ }}}}{T_1}} \hfill \cr {{S_2} = {1 \over 2}{A_1}T_2^2 + {F_1}{T_2}} \hfill \cr {{S_3} = - {1 \over 6}{J_{max}}T_3^3 + {1 \over 2}{A_2}T_3^2 + {F_2}{T_3}} \hfill \cr {{S_5} = - {1 \over 6}{J_{max}}T_5^3 + {V_{max}}{T_5}} \hfill \cr {{S_6} = - {1 \over 2}{A_5}T_6^2 + {V_5}{T_6}} \hfill \cr {{S_7} = {1 \over 6}{J_{max}}T_7^3 - {1 \over 2}{A_6}T_7^2 + {V_6}{T_7}} \hfill \cr {{S_4} = {S_{{\rm{seg }}}} - \left({{S_1} + {S_2} + {S_3} + {S_5} + {S_6} + {S_7}} \right)} \hfill \cr } } \right. where: S_1-7 – displacement during each phase, T_1-7 – duration of each phase, which is computed as: 20 ${\begin{array}{l} T_{1} = T_{3} = T_{5} = T_{7} = \frac{A_{m a x}}{J_{m a x}} \\ T_{2} = \frac{V_{m a x} - V_{s t a r t}}{A_{m a x}} - T_{1} \\ T_{6} = \frac{V_{5} - V_{6}}{A_{m a x}} \end{array}$ \left\{ {\matrix{ {{T_1} = {T_3} = {T_5} = {T_7} = {{{A_{max}}} \over {{J_{max}}}}} \hfill \cr {{T_2} = {{{V_{max}} - {V_{start }}} \over {{A_{max}}}} - {T_1}} \hfill \cr {{T_6} = {{{V_5} - {V_6}} \over {{A_{max}}}}} \hfill \cr } } \right.

To plan the S-Curve feedrate profile for a NURBS toolpath, the critical points of the NURBS curve have to be identified. These points mark fragments of the path with curvature maxima. At these points, the feedrate profile should achieve a minimum to improve interpolation and tracking accuracy. To select the critical points, a chord error feedrate limit function is defined for each sampling point[12]: 21 $v_{c h o r d} (u_{i}) = m a x (\frac{2}{Δ t} \sqrt{ρ^{2} (u_{i}) - {(ρ (u_{i}) - ϵ_{m a x})}^{2}}, V_{m a x})$ {v_{chord{\rm{ }}}}\left({{u_i}} \right) = max\left({{2 \over {\Delta t}}\sqrt {{\rho ^2}\left({{u_i}} \right) - {{\left({\rho \left({{u_i}} \right) - {_{max}}} \right)}^2}} ,{V_{max}}} \right) 22 $κ (u_{i}) = \frac{1}{ρ (u_{i})} = \frac{‖ C^{'} (u_{i}) \times C^{''} (u_{i}) ‖}{{‖ C^{'} (u_{i}) ‖}^{3}}$ \kappa \left({{u_i}} \right) = {1 \over {\rho \left({{u_i}} \right)}} = {{\left\| {{C^\prime }\left({{u_i}} \right) \times {C^{\prime \prime }}\left({{u_i}} \right)} \right\|} \over {{{\left\| {{C^\prime }\left({{u_i}} \right)} \right\|}^3}}} where: κ(u_i), ρ(u_i) – curvature and curvature radius for parameter u_i at sample point i, ϵ_max – maximum allowable chord error, V_max – maximum allowable feedrate, i – index of current feedrate sampling point.

Chord error is the distance between a line connecting consecutive interpolation points and the actual curve approximated locally by a circle with a radius equal to the local curvature radius of the curve. Chord error is the measure of interpolation accuracy and v_chord(u_i) determines the maximum feedrate at each sampling point between start and end of the curve, for which the chord error limit is not violated. Minima of the chord error feedrate limit function determine sections of the NURBS toolpath between which S-Curve feedrate profiles are generated. Values of the v_chord function determine initial and final velocities of the S-curve segment.

The S-curve segments are planned for each segment with maximum predefined feedrate, tangential acceleration and jerk so that the distance covered by each profile segment (S_seg) is equal to the arc-length distance between critical points. If the distance between the critical points is too short to fit the full s-curve, the maximum feedrate in that segment is iteratively decreased until the profile displacement fits the arc-length distance between the critical points.

3.2.

Feedrate profile optimization using PSO

Feedrate optimization is performed on a feedrate profile in the form of a B-Spline curve. The B-spline is similar to the NURBS curve but without weights. The B-Spline is defined as follows: 23 $C (u) = \sum_{i = 1}^{k} N_{i, n} (u) P_{i}$ {\bf{C}}({\bf{u}}) = \sum\nolimits_{{\rm{i}} = 1}^{\rm{k}} {{{\rm{N}}_{{\rm{i}},{\rm{n}}}}} ({\rm{u}}){{\bf{P}}_{\bf{i}}}

The S-curve polynomial feedrate profile, generated in the previous step, is converted to a one-dimensional cubic B-Spline representation. The feedrate profile is reparametrized from a function of time to a function of the NURBS toolpath parameter.

This procedure is performed in the following steps:

–
For each segment delimited by the critical points of the NURBS curve the initially planned S-curve feedrate profile is sampled and the with equal time intervals. Ten samples per segment were used with the number being determined by trial and error to balance accuracy and computational complexity.
–
The toolpath is interpolated and for each sampled feedrate profile value the corresponding value of NURBS curve parameter u is also sampled.
–
Using the curve parameter and feedrate samples recorded in the previous step, a cubic B-spline curve is fitted in that segment using least-squares. The mean squared error between sampled feedrate points and value of the B-Spline evaluated for corresponding parameter values are minimised [6]. The endpoints of the fitted segment are equal to the critical points at the start and end of the segment.
–
The above steps are repeated for each segment and the final feedrate profile consists of a joined B-spline parametrized with the NURBS toolpath parameter.
–
The reparametrized B-Spline profile is returned as a series of control points

The optimization of cubic B-Spline feedrate profile is performed by adjusting the control point values, which are treated as optimization variables.

The objective function, can be formulated as: 24 $f (P) = \frac{1}{N} \sum_{i = 1}^{N} (1 - \frac{P_{i}}{V_{m a x}})$ {\bf{f}}({\bf{P}}) = {1 \over N}\sum\nolimits_{{\rm{i}} = 1}^{\rm{N}} {\left({1 - {{{{\bf{P}}_{\bf{i}}}} \over {{V_{max}}}}} \right)} where: N – number of optimized control points, P_i – control point, i – control point index, V_max – maximum allowable feedrate.

The objective function maximizes the normalized sum of control points so that the feedrate is as high as possible.

The optimization problem is subject to constraints of maximum axial velocity, acceleration and jerk: 25 $c_{v (1, 2)} = \frac{| v_{(1, 2)} (u_{i}) |}{v_{m a x}} - 1$ {c_{v(1,2)}} = {{\left| {{v_{(1,2)}}\left({{u_i}} \right)} \right|} \over {{v_{max}}}} - 1 26 $c_{a (1, 2)} = \frac{| a_{(1, 2)} (u_{i}) |}{a_{m a x}} - 1$ {c_{a(1,2)}} = {{\left| {{a_{(1,2)}}\left({{u_i}} \right)} \right|} \over {{a_{max}}}} - 1 27 $c_{j (1, 2)} = \frac{| j_{(1, 2)} (u_{i}) |}{j_{m a x}} - 1$ {c_{j(1,2)}} = {{\left| {{j_{(1,2)}}\left({{u_i}} \right)} \right|} \over {{j_{max}}}} - 1

Where: c_{v(1, 2)}, c_{a(1, 2)}, c_{j(1, 2)} – normalized constraint violation values for velocity, acceleration and jerk in each axis (1 or 2). These constraints are computed based on axial velocity, acceleration and jerk computed using equations (16–18). v_max, a_max, j_max – maximum allowable value of axial velocity, acceleration and jerk.

Minimisation of the objective function is performed using the widely proven Particle Swarm Optimization algorithm (PSO)[14]. PSO is a swarm intelligence, nature-inspired optimization method in which virtual particles (candidate solutions) move through the search space and converge on the minimum of the objective function, which is indicated by the “global best” position (last position with the best fitness so far). In each iteration, the velocities and positions of each particle (solutions) are updated according to the following equations [13]: 28 $v_{i, j}^{*} = ω \cdot v_{i, j} + ϕ_{1} \cdot r a n d_{(0, 1)} \cdot (p_{i, j} - x_{i, j}) + ϕ_{2} \cdot r a n d_{(0, 1)} \cdot (g_{i} - x_{i, j})$ v_{i,j}^* = \omega \cdot{v_{i,j}} + {\phi _1}\;\cdot\;r\;an{d_{(0,1)}}\;\cdot\;\left({{p_{i,j}} - {x_{i,j}}} \right) + {\phi _2}\;\cdot\;r\;an{d_{(0,1)}}\;\cdot\;\left({{g_i} - {x_{i,j}}} \right) 29 $x_{i, j}^{*} = x_{i, j} + v_{i, j}^{*}$ x_{i,j}^* = {x_{i,j}} + v_{i,j}^*

Where: x_{i, j} – position of particle j in step i, v_{i, j} – position of particle j in step i, P_{i, j} – personal best position of particle j in step i, g_i - global best position at iteration i, ω, φ₁, φ₂, – parameters of the algorithm (inertia weight, cognitive coefficient, social coefficient), rand_{(0, 1)} – random number between 0 and 1 drawn from a uniform distribution.

The standard Particle Swarm Optimization method does not handle arbitrary constraints. Only the optimization variables are limited. To handle axial velocity, acceleration and jerk constraints, such defined by equations (23–25), the method has to be extended with some constraint handling technique. There are several constraint handling techniques such as Deb’s method, Penalty method and Augmented Lagrangian Method [18].

The Augmented Lagrangian Method, used in this paper, adds a penalty factor to the objective function which is dependent on the amount of constraint violation. This augmented objective function creates the Augmented Lagrangian which is defined as [15]: 30 $L (P) = f (P) + \frac{ρ}{2} {\sum_{i = 1}^{M} (\max {0, c_{i} (P) + \frac{λ_{i}}{ρ}})}^{2}$ L({\bf{P}}) = f({\bf{P}}) + {{\rm{\rho }} \over 2}{\sum\nolimits_{i = 1}^M {\left({\max \left\{ {0,{c_i}({\bf{P}}) + {{{{\rm{\lambda }}_i}} \over {\rm{\rho }}}} \right\}} \right)} ^2} where: f(P)-objective function value for current solution vector P, ρ – penalty factor, c_i(P) – value of constraint i for current solution vector P, λ_i – Lagrange multiplier for constraint i, M - number of constraints.

The Augmented Lagrangian is minimised by the Particle Swarm Optimization several times. After each minimisation, the Lagrange multipliers and penalty factor are updated to find their optimal value. This leads to the solution of the unconstrained problem that minimises the augmented Lagrangian to be equivalent to the solution of the original constrained problem [17].

The update rules of the penalty factor and Lagrange multipliers are defined as [15]: 31 $λ_{i}^{*} = m i n {\max {0, λ_{i} + ρ c_{i} (P)}, 10^{20}}$ {\rm{\lambda }}_i^* = min\left\{ {\max \left\{ {0,{{\rm{\lambda }}_i} + {\rm{\rho }}{c_i}({\bf{P}})} \right\},{{10}^{20}}} \right\} 32 $ρ_{i}^{*} = {\begin{matrix} 10 ρ & I C M^{*} > 0.5 * I C M \\ ρ & I C M^{*} \leq 0.5 * I C M \end{matrix}$ {\rm{\rho }}_i^* = \left\{ {\matrix{ {10\rho } & {IC{M^*} > 0.5*ICM} \cr \rho & {IC{M^*} \le 0.5*ICM} \cr } } \right. 33 $IC M * = \max {I C M, | \max {c_{i} (P), - \frac{λ_{i}}{ρ}} |}$ {\rm{IC}}M* = \max \left\{ {ICM,\left| {\max \left\{ {{c_i}({\bf{P}}), - {{{\lambda _i}} \over \rho }} \right\}} \right|} \right\}

The Infeasibility Compatibility Measure (ICM) is the measure of constraint violation. The idea behind the update method is to increase the penalty factor and Lagrange multipliers if the ICM increases which should lead to the elimination of the constraint violation in subsequent runs of the optimization. The result of each augmented Lagrangian minimisation is used as an initial guess for the next cycle.

Optimization of the B-Spline feedrate profile is performed in a moving window approach. A constant number of control points is optimized during each run, equal to 50. Each window has an overlap o 10 control points with the previous window. The first and last control point in the window is not modified. In each optimization window the constraints are evaluated at equally spaced constraint evaluation points (CEP) with the distance set to 0.5mm.

Initialization of the Particle Swarm Optimization algorithm is usually done by randomizing the swarm in the whole search space. In this approach the pre-generated feedrate profile is used to initialize one particle to exactly match the initial solutions. Other particle positions are randomized around that initial solution. This narrows down the search space and improves convergence speed.

The algorithm parameters suggested in [13], commonly used in literature and in this paper are: ω = 0,7298, ϕ₁ = ϕ₂ = 1,4926. Other important parameters include the number of particles and termination conditions. The termination conditions are usually the number of iterations or function evaluations. The number of particles was chosen to be 25 and the termination criteria for a single PSO optimization is chosen as 200 iterations. The number of maximum augmented Lagrangian optimizations is chosen as 25.

4.

EXPERIMENTAL VERIFICATION OF THE FEEDRATE GENERATION ALGORITHM

4.1.

Experimental setup

The experimental setup, used to experimentally verify the proposed feedrate optimization algorithm, consists of an H-Bot gantry numerically controlled machine, shown in (Fig. 6), with a PC-based numerical controller. The CNC machine is driven by Delta ASDA-A3-E 750 W servo drives with permanent magnet servo motors. The HTD toothed belt is actuated by pulleys driven by the servo motors with 1:5 reduction gearboxes. The drives are controlled by the PC-based controller using the EtherCAT communication bus.

The PC based controller is implemented on an industrial PC with the Linux operating system with RT-Preempt real-time kernel. An open-source IGH EtherCAT Master stack is implemented as a real-time module to communicate with the servo drives. The communication is conducted according to the CAN in Automation (CiA) 402 standard drive profile. The real-time thread cycle that governs all real-time models, including the communication cycle, was set to 250μs.

The feedrate optimization algorithm is implemented as a non-real-time user program. The generated optimized feedrate profile is then sent to the real-time NURBS interpolator. The interpolator generates the consecutive position commands to the drives in real-time during machine operation. The data such as actual position, speed, current etc. is collected from the drives and linear encoders in real-time and saved to file through a FIFO buffer. The data can be later analysed in tools such as Matlab. The experimental setup is described in more detail in [16]. A block schematic of the experimental setup’s control system is presented in (Fig. 7).

4.2.

Experimental results

In order to experimentally verify the proposed method an example NURBS toolpath (“flower”) was used (Fig. 2) which has a total length of 7001.5 mm. First an initial guess s-curve was generated with limits of feedrate, tangent acceleration and tangent jerk were set to: 1250 mm/s, 12000 mm/s2 and 120000 mm/s3 respectively. The generated initial S-curve profile is shown in (Fig. 8).

When generating the initial profile, the axial velocity, acceleration and jerk limits are not taken into account directly, but can be limited indirectly by limiting feedrate, tangent acceleration and jerk. The maximum feedrate was set so that the axial velocity limits are within the nominal motor velocity (314 rad/s or 3000 rev/min). Expected axial acceleration and jerk limits are 1700 rad/s2 and 51500 rad/s3. The acceleration limit serves to limit drive current to its nominal value of 5A. Jerk limit serves to limit vibration and is set experimentally. Acceleration and jerk limits are occasionally violated, which could result in motor overloads, following error spikes or excessive vibration. The generated trajectory run time is 7.821s.

The generated initial feedrate profile was converted to a cubic B-Spline and input to the feedrate optimization algorithm. The obtained profile had 402 control points. The optimization algorithm was then run on the B-spline profile. The trajectory execution time is 7.667s while the average optimization time from 10 separate runs is 76.1s. The resulting trajectory was executed on the H-Bot CNC machine using the “flower” toolpath. The toolpath with indicated feedrate is shown in (Fig. 9). The optimal feedrate profile with axial velocity, acceleration and jerk for both axes is shown in (Fig.10). Measured drive current and following errors for both axes are shown in (Fig.11).

It can be seen that the generated profile does not violate axial constraints despite the feedrate profile reaching higher values than the initial 1250mm/s. The resulting trajectory also has a shorter runtime than the initial trajectory. The axial acceleration and jerk could be limited indirectly by limiting tangent acceleration and jerk, but this would result in sub-optimal execution times and lower productivity. Replacing the standard S-Curve with a cubic B-Spline feedrate profile caused the jerk to be continuous, which improves motion smoothness and reduces potential vibration. The drive current is within nominal values and the following error is bounded.

5.

CONCLUSION

This paper presents a feedrate optimization algorithm for a CNC machine with non-Cartesian H-Bot kinematics and Non-Uniform Rational B-Spline (NURBS) toolpaths. The feedrate profile is optimized using the Particle Swarm Optimization (PSO) algorithm. To speed up convergence and improve reliability, the feedrate profile is first planned as an S-Curve and then converted to a B-Spline curve. The B-Spline profile is then adjusted by the PSO algorithm with the Augmented Lagrangian constraint handling technique to maximize feedrate.

Experimental results show that the velocity, acceleration and jerk in the machine axes are within their prescribed limits. Compared to the initial profile, the trajectory execution time was reduced, and no constraint violations were observed. Axial velocity, acceleration and jerk limitation also limit the drive current to nominal values and can serve to indirectly limit the following error.

Future research will focus on directly handling following and contour error constraints and extension of the proposed method to other non-Cartesian kinematics.

Jerk Limited Feedrate Profile Optimization for Nurbs Toolpaths in Cnc Machines with H-Bot Kinematics

Full Article

Paradigm

My account