Toward Optimal Robot Machining considering the Workpiece Surface Geometry in a Task-Oriented Approach

1. Introduction

Nowadays, robotics is an indispensable technology in many industries, especially in manufacturing, since it represents a major building block for fully automated production lines such as in the automotive industry. Typically, standard applications of industrial robots are designed such that they perform well-defined repetitive tasks of manipulation, assembly, palletizing, welding, painting, etc [1]. However, on the other hand, a shift from high volume/low mix to low volume/high mix and custom based production impacts manufacturing seriously. It requires robots which are easier to install, program and operate, in order to increase robotization in small and medium-sized enterprises, where collaborative and more flexible automation appears to be a useful option. Thus, collaborative robots (cobots), not only provide safe physical coexistence and interaction and cooperation in a common workspace with human workers [2,3], but are also easier to use, and are gaining popularity in industries with a growing demand for highly customized products, since they increase manufacturing flexibility, filling the gap between fully automated systems and manual production [4,5]. However, in order to guarantee safe physical human-robot interaction, collaborative robots are designed with substantially reduced capabilities of speed and force/torque in comparison with standard industrial robot arms [6].

Beside the above mentioned industrial robot applications there are also others, such as more complex robot machining [7] or milling [8], deburring [9], surface finishing [10], grinding [11], polishing [12,13], hammer-peening [14,15], etc. In some processes the technology requires a significant machining force in the tool feed direction, and in some the machining force is required in the direction only toward the workpiece surface. The technological process of metal surface treatment with hammer-peening, which is applied in the tool- or mold- or die-making industry, is such a cold forging process, in which a ball made of a carbide solid is struck with a high frequency on the metal surface of the workpiece by a micro-forging hammer tool [16]. Hammer-peening improves the smoothness of the surface, hardens it, and eliminates internal stresses. In this manufacturing process, the machine hammer-peening tool (MHP tool) is moving freely along the surface of the workpiece (that may have a complex geometry) following the machining path, which connects all surface points of the patch area to be forged. The working angle of the MHP tool to the workpiece is usually determined as perpendicular to its surface along the direction of travel of the tool. However, such surface finishing process is complex and thus a very wasteful process, which has a large impact on the cost of the overall processing in the tool-making industry. It demands to combine both a collaborative, intelligence-based and a cooperative human-robot-based technological approach, as in the robotic polishing application [17]. Thus, it may involve not only robotized solutions, e.g., such as [18,19], but also the collaborative automation with a proper flexibility due to the customized high mix/low volume nature of the production type in the tooling industry. Therefore, cobots can be considered as a justified replacement for traditional industrial robots in such applications.

However, due to low power built-in cobot actuators (in comparison with standard industrial robot arms) their introduction in a machining process requires more careful trajectory planning, in order to assure the feasibility of the robot task, especially in the case of manufacturing processes with complex continuous paths where the complexity of robot path planning increases significantly [20]. Optimal relative workpiece/robot placement and robot path/trajectory planning considering this issue thus become even more important, in order to provide rapid set-up of a robotic system in flexible high mix/low volume applications. Such planning problems have been interesting research topics for many years, and their application in practice still attracts attention in the research community. The researchers optimized the location of the robot for generating maximum task-space velocity with minimum joint velocities [21], and robot-to-workpiece placement for large scale welding systems [22], generated a kinematic performance map based on a kinetostatic condition index that was used to optimize robot configurations in a polishing application [23], introduced a custom index for robot base placement optimization demonstrated in a trim application in shoe manufacturing [24], and optimized a workpiece placement for the robotic operation in challenging manufacturing tasks [25,26] and surface finishing [20,27]. An interesting new optimization approach was also introduced, to maximize the available velocities of the end-effector during a task execution of path following in robot machining called the Decomposed Twist Feasibility method [28].

There have been several attempts at trajectory optimization for an industrial robot with different objectives. However, the main well-known basic concept that enables analysis of robot motion capability is called manipulability [29]. The pioneering work was presented by Yoshikawa in 1985 [30], who introduced a qualitative method for assessment of robot motion capabilities based on the so-called manipulability ellipsoid, that was derived from a robot velocity kinematics description by the application of the singular value decomposition of the associated Jacobian matrix and Euclidean norm metrics. The introduced manipulability index was proportional to the ellipsoid volume, and it should measure how easy or difficult it is for a robot to move in its Cartesian operational workspace, and it may also represent a distance to the singular configuration of the robotic arm. Beside the volume of the ellipsoid, other indices also derived from the ellipsoid appeared later, such as minimal singular value and condition number [31]. However, the operational space twist vector involves different units of linear and angular velocity; beside this, robot joints can also be of different types. These complicate the manipulability analysis from the point of view of physical consistency [32,33]. Thus, in order to avoid the problem of the dimensional dependence when both position and orientation of the robot end-effector are included in kinematic equations, new manipulability indices were introduced, which were based on the introduction of auxiliary points on the end-effector that provided additional linear velocity information and a redundant formulation of the velocity equations, instead of the combination of linear and angular velocities [34]. However, in this point-based approach the determination of the auxiliary points is arbitrary, and the measure of the obtained dimensionally homogeneous Jacobian matrix is not invariant with respect to changes of the auxiliary points used to express the end-effector velocity. A general approach to the problem of dimensional nonhomogeneous matrices in the velocity kinematics description presents the introduction of weighted norms [35], and the selection of the corresponding weighting matrices to set translational components of the Jacobian matrix in relation with the rotational components [36] is again arbitrary, e.g., manual selection of the weights in a task dependent measure [37], or their selection can be based on the computation of the minimal principal angle between the translational and rotational subspace in a task-oriented approach [38]. However, it has been shown that all such manipulability indices suffer from the nonexistence of “natural” metrics, and therefore from non-invariance in the sense of the choice of the selected artificial metric functions, which is arbitrary, employed in their definition [39,40]. Although the translational and rotational operational velocities can be separated by exploiting and extending of the concept of manipulability in both weak and strong senses [41], the problem remains. Another possibility to overcome the problem of a non-homogeneous Jacobian matrix is based on the apparent power concept of the robot mechanism, resulting in a homogeneous formulation of the problem, regardless of having mixed units in the velocity kinematics description [42]. However, the relation of the proposed measure to the relevant robot operational quantity such as velocity is unclear. An alternative to manipulability ellipsoids is presented by manipulability polytopes [28,43,44,45,46,47,48], which can provide more accurate information about the operational space motion capability, since they consider infinity norm metrics instead of Euclidean norm metrics. However, the problem of dimensional a non-homogeneous Jacobian matrix due to the unit inconsistency in robotics still present a challenge, since most approaches demonstrated limitations in terms of their physical interpretations, though task-oriented homogeneous Jacobians and associated performance indices showed some promising results and potential for further development [28,49].

All the performance indices above and many others not mentioned are based on the Jacobian matrix of the velocity kinematics description. In this paper, we focus on the problem of path following velocity performance in robot machining of a workpiece with complex geometry, similarly as in our previous work [28] , where we considered a workpiece with a predefined path. Now, we propose a task-oriented robot kinematics description, which considers motion constraints imposed by the workpiece surface geometry explicitly; thus, a priori knowledge about the machining path is not required to evaluate the available velocity performance. We show that the associated manipulability ellipsoid is reduced in dimension, such that we can perform velocity planning solely within the translational operational space projected on the surface tangent plane, and thus the problem is avoided of the nonhomogeneous twist space due to the mixed units. The proposed constrained kinematics description can provide a basis, not only for workpiece/robot placement optimization, but also for optimal path planning based on kinematic manipulability. We demonstrate the proposed approach by a numerical simulation case study on two different workpiece examples.

The rest of the paper is organized as follows: in Section 2 we provide the necessary background information from the surface differential geometry, in Section 3 we develop the proposed approach, in Section 4 we show the numerical experiments, and section 5 concludes the paper.

2. Background

In this Section we provide some basic information about the differential geometry of surfaces and curves on surfaces which are relevant to this paper.

Let us consider surface $\mathrm{S}$ in a three-dimensional Euclidean space that is parametrized by the position vector $r$

$$r=\left[x(u,v),y(u,v),z(u,v)\right],$$

(1)

as a function in the parametric form

$$r=r(u,v),$$

(2)

where $u$ and $v$ are independent parameters in a closed rectangle [50,51]. In this paper, we consider an explicit surface such that the z-coordinate can be expressed as a function of both the x- and y-coordinates:

$$z=f(x,y),$$

(3)

where $f$ is at least twice the differentiable real-valued function, and the parametric form is derived by setting $x=u$ and $y=v$, with $\left[x,y\right]\in \mathit{D}$, where $\mathit{D}$ is a bounded connected domain in the real xy-plane. The surface is then given by

$$r=\left[x,y,f(x,y)\right].$$

(4)

A curve on the surface is given by the parametrization $x=x(t)$, $y=y(t)$, and $z=f(x(t),y(t))$, where $t\in \left[0,T\right]$, and with $\left[x,y\right]\in \mathit{D}$. Then $r=r(t)$ is a curve lying on and embedded in the surface (4). The tangent vector to the curve on the surface is evaluated by differentiating $r(t)$ with respect to the parameter $t$ using the chain rule, and is given by:

$$\dot{r}(t)=r_x\dot{x}+r_y\dot{y},$$

(5)

where the subscripts x and y denote partial differentiation with respect to x and y such that $r_x=\partial r/\partial x$ and $r_y=\partial r/\partial y$, respectively, and the dot denotes differentiation with respect to the parameter $t$, such that $\dot{x}=dx/dt$ and $\dot{y}=dy/dt$. Note that (5) can be written in a form independent of the choice of parameter:

$$dr=r_{x}dx+r_{y}dy,$$

(6)

where $d(.)$ denotes a differential. The differential arc length of the curve $ds$ is given as:

$$ds=\left|\frac{dr}{dt}\right|dt=\left|\dot{r}\right|dt=\sqrt{\dot{r}\circ \dot{r}}dt,$$

(7)

where $\circ$ stands for the vector dot product operator.

Let $P=P(x,y)$ be a point on the regular surface $\mathit{S}$. Then $r_x$ and $r_y$ are two independent surface tangent vectors at point $P$, which span a tangent plane. The tangent plane at point $P$ on surface $\mathit{S}$ can be considered as a union of all tangent vectors which can be formed as a linear combination of $r_x$ and $r_y$. The unit normal vector $\widehat{n}$ on the surface at point $P$ can be defined as

$$\widehat{n}=\widehat{n}(x,y)=\frac{r_x\times r_y}{‖r_x\times r_y‖}.$$

(8)

The unit normal vector is mapped from the tangent vectors $r_x$, $r_y$. It is perpendicular to the tangent plane, and, obviously, it is also orthogonal to both $r_x$ and $r_y$. Thus, at any point on the surface we have a triple of vectors $r_x$, $r_y$, and $\widehat{n}$, that are linearly independent on the regular surface $\mathit{S}$. Note that the differential change of the unit normal vector can be written as:

$$d\widehat{n}={\widehat{n}}_{x}dx+{\widehat{n}}_{y}dy,$$

(9)

where ${\widehat{n}}_x=\partial \widehat{n}/\partial x$ and ${\widehat{n}}_y=\partial \widehat{n}/\partial y$ denote partial differentials of the unit normal vector with respect to x and y, respectively.

The first fundamental form describes the way of measuring the distances on a surface, i.e., a surface metric. It determines the arc length of a curve on the surface, and is defined as [51]:

$${\rm I}:d{s}^2=dr\circ dr.$$

(10)

If we consider (6) then it can be derived as follows

$${\rm I}(dx,dy):d{s}^2=Ed{x}^2+2Fdxdy+Gd{y}^2,$$

(11)

where the coefficients $E=r_x\circ r_x$, $F=r_x\circ r_y$, and $G=r_y\circ r_y$ are called the coefficients of the first fundamental form. Due to its quadratic-bilinear form of the coordinates’ differentials on the surface, it is often presented by the symmetric matrix:

$${\rm I}=\left[\begin{array}{cc}E& F\\ F& G\end{array}\right],$$

(12)

which is a positive definite, i.e., $EG-F^2>0$ at regular points on the surface.

The second fundamental form characterizes the local structure of the surface shape in a neighborhood of a point. It describes how the surface deviates from the tangent plane. It can be defined as [51]

$${\rm I}{\rm I}:(-dr\circ d\widehat{n}),$$

(13)

It can be derived as

$${\rm I}{\rm I}(dx,dy):Ld{x}^2+2Mdxdy+Nd{y}^2,$$

(14)

where the coefficients $L=-r_x\circ {\widehat{n}}_x=r_{xx}\circ \widehat{n}$, $M=-r_x\circ {\widehat{n}}_y=-r_y\circ {\widehat{n}}_x=r_{xy}\circ \widehat{n}=r_{yx}\circ \widehat{n}$, and $N=-r_y\circ {\widehat{n}}_y=r_{yy}\circ \widehat{n}$ are called the coefficients of the second fundamental form. The matrix of the second fundamental form can be read as:

$${\rm I}\text{}{\rm I}=\left[\begin{array}{cc}L& M\\ M& N\end{array}\right].$$

(15)

The differential $d\widehat{n}$ is called a Weingarten map, that describes the change in the normal direction as we move from one point to another. It can be expressed in terms of the first derivatives of the position vector, and the partial derivatives of the unit normal vector can then be expressed in terms of the basis $\left\{r_x,r_y\right\}$ [51,52]:

$$\left[\begin{array}{cc}{\widehat{n}}_x& {\widehat{n}}_y\end{array}\right]=-\left[\begin{array}{cc}{r}_x& r_y\end{array}\right]W,$$

(16)

where $W$ is the matrix, that can be calculated by (12) and (15) as:

$$W={{\rm I}}^{-1}\cdot {\rm I}\text{}{\rm I}=\frac{1}{EG-F^2}\left[\begin{array}{cc}GL-FM& GM-FN\\ EM-FL& EN-FM\end{array}\right],$$

(17)

where the operator $\cdot$ stands for matrix multiplication. The result is known as Weingarten equations:

$${\widehat{n}}_x=\frac{FM-GL}{EG-F^2}{r}_x+\frac{FL-EM}{EG-F^2}{r}_y,$$

(18)

$${\widehat{n}}_y=\frac{FN-GM}{EG-F^2}{r}_x+\frac{FM-EN}{EG-F^2}{r}_y.$$

(19)

which map from the tangent plane to the tangent plane, i.e., the vectors ${\widehat{n}}_x$ and ${\widehat{n}}_y$ are expressed as a linear combination of the tangent vectors $r_x$, $r_y$, respectively. Thus, they lie in the tangent plane, and are orthogonal to the unit normal vector $\widehat{n}$.

An important surface characteristic is also its curvature, which shows us how much the surface bends or deviates from a flat surface. The curvature measure should involve the rate at which surface $\mathit{S}$ leaves the tangent plane to $\mathit{S}$ at point $P$. A crucial tool to measure the curvature on a surface is the Weingarten map, which contains complete information about a surface’s curvature. It can be used to derive scalar measures to show how ‘curved’ the surface is at some point. There are several measures for the curvature of a surface in ${ℝ}^3$. Gauss curvature $K$ and Mean curvature $H$ on surface $\mathit{S}$ at point $P$ can be formulated by (20) and (21), respectively.

$$K=\det (W)=\frac{LN-M^2}{EG-F^2}$$

(20)

$$H=\frac{1}{2}\mathrm{trace}(W)=\frac{EN+GL-2FM}{2(EG-F^2)}$$

(21)

Other measures, such as normal curvature and principal curvatures, are also possible [52].

3. Methodology

3.1. Preliminaries

The serial robot arm manipulator can be described as an open chain mechanism connecting links and joints having n- degrees-of-freedom (DOF). The base link is normally considered fixed to the ground, and the final link, with mounted end-effector, may operate freely in the robot task space. To provide all possible positions and orientations within its 3D workspace, the manipulator must have proper mechanism geometry with at least 6 DOFs. The velocity mapping from the joint space to the operational space can be described as the forward velocity kinematics model [29]:

$$\nu =J(\theta )\dot{\theta },$$

(22)

where $J\in {ℝ}^{6xn}$ is the geometric Jacobian matrix in the robot`s base frame, $\theta \in {ℝ}^n$ and $\dot{\theta }\in {ℝ}^n$ are the vectors of joint positions, which define the current mechanism configuration, and their velocities, respectively. Furthermore, we denote the so-called twist vector with $\nu \in {ℝ}^6$

$$\nu =\left[\begin{array}{c}v\\ \omega \end{array}\right],$$

(23)

which represents the combination of two physically different components: $v\in {ℝ}^3$ is the vector of the end-effector linear velocity, and $\omega \in {ℝ}^3$ is the vector of the end-effector angular velocity. The inverse velocity kinematics problem can be defined as finding a proper solution of the joint velocity $\dot{\theta }$ for the given twist $\nu$. Note that the twist vector is not homogeneous in physical units.

In the following, we consider a nonredundant manipulator with $n=6$ in a non-singular configuration. Although in general the joints can be of different types, we assume rotational joints without a significant loss of generality, since most modern robot manipulators use only rotary joints with one DOF.

The most widely used kinematic measure that refers to manipulator manipulability and indicates how dexterous a robot is at a given configuration $\theta$ is derived from the condition [30]:

$${\dot{\theta }}^T\dot{\theta }\le 1,$$

(24)

which maps the unit sphere in the joint velocity space to the velocity manipulability ellipsoid in the operational task space:

$${\nu }^T{(J{J}^T)}^{-1}\nu \le 1.$$

(25)

Here, we consider that the solution to the inverse velocity kinematics can be expressed by the full rank Jacobian as:

$$\dot{\theta }=J^{-1}\nu .$$

(26)

We can factorize the Jacobian matrix by Singular Value Decomposition (SVD) into

$$J=U\Sigma V^T,$$

(27)

where $\Sigma$ is a diagonal matrix with positive singular values ${\sigma }_1\ge {\sigma }_2\ge \mathrm{...}\ge {\sigma }_6>0$ on the diagonal, and $U$and $V$ are unitary matrices representing the orthonormal base of the Jacobian column space $Col(J)$and the orthonormal base of the Jacobian row space $Row(J)$, respectively [53,54]. If we insert (27) into (25), it yields:

$${(U^T\nu )}^T{\Sigma }^{-2}(U^T\nu )\le 1,$$

(28)

where ${\Sigma }^{-2}$ is a diagonal matrix with ${\sigma }_1^{-2},{\sigma }_2^{-2},\mathrm{...},{\sigma }_6^{-2}$ on the diagonal. One can introduce a new variable $w=U^T\nu$, $w={\left[w_1,w_2,\mathrm{...},w_6\right]}^T$; then (28) can be rewritten as:

$$w^T{\Sigma }^{-2}w={\displaystyle \sum _{i=1}^6{(\frac{w_i}{{\sigma }_i})}^2}\le 1,$$

(29)

which describes an ellipsoid geometrically with the axes` directions defined by the columns of the unitary matrix $U$, whereas the singular values ${\sigma }_i$ determine the axes` lengths. The manipulability ellipsoid shows how far away the manipulator is from a singular configuration in a certain direction in the task space: if the manipulator is near the singular configuration, then the dexterity is low, and contrarily, if the manipulator is far from the singular configuration, then the dexterity is better. Using the manipulability ellipsoid one can quantify how close a given robot posture is to a singularity. Note that the ellipsoid is described in the 6-dimensional Euclidian space, which is not suitable for visual geometrical representation, and is not intuitive.

Although the manipulability ellipsoid is quite popular among the research community since it enables the derivation of numerous manipulability indices [55], it suffers from serious limitations, that arise due to the physical non-homogeneity of the twist vector and physical inconsistency of the Jacobian matrix [32,41].

3.2. The end-effector motion constraints

During the workpiece machining the tool tip must follow a certain path embedded to the workpiece surface. The tool path is constrained by position and orientation. In our case, we assumed that the tool must maintain orthogonal orientation to the curved workpiece surface, i.e., the tool orientation must be adjusted continuously while traveling along the path. In this Section, we derive velocity synchronization, which arises due to the constraints described above.

The tool tip must follow the position of the path embedded on the workpiece surface. Let a particle refer to the tool tip. If a particle is moving in time along the path constrained on the workpiece surface, then the particle velocity can be expressed as:

$$v=\dot{r},$$

(30)

If we consider (5), then the tool tip velocity can be derived as

$$v=\left[\begin{array}{cc}{r}_x& r_y\end{array}\right]\left[\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right].$$

(31)

Note that the tool linear velocity vector $v$ is a linear combination of the tangent vectors $r_x$ and $r_y$. Thus, it lies in a tangent plane of the workpiece surface. On the other side, (31) it can be read in a way that $v$ is to be synchronized with the velocity vector in the xy-plane ${\left[\begin{array}{cc}\dot{x}& \dot{y}\end{array}\right]}^T$ by the matrix $S_v$:

$$S_v=\left[\begin{array}{cc}{r}_x& r_y\end{array}\right],$$

(32)

such that we can express:

$$v=S_v\left[\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right],$$

(33)

where $S_v\in {ℝ}^{3\times 2}$ is with full rank, i.e. $\mathrm{rank}(S_v)=2$.

The tool orientation must be aligned by the normal vector perpendicular to the curved workpiece surface. Let a particle refer to the tool tip. If a particle is traveling along the path constrained on the workpiece surface, then, in each point of the path, the normal vector associated to the particle position can have a different orientation. Not only that the particle is translating, but also that the associated surface normal vector can rotate synchronously in this case. We can observe the velocity of the unit normal vector $\widehat{n}=\widehat{n}(t)$ that can be described as:

$$\dot{\widehat{n}}=\omega \times \widehat{n},$$

(34)

where $\times$ refers to the vector cross product operator. If we assume reasonably that the normal vector rotation about itself is zero, then the rotation axis is orthogonal to the unit normal vector, and it lies in the surface tangent plane of the observed point on the path. Then the angular velocity of the unit normal vector can be derived as:

$$\omega =\widehat{n}\times \dot{\widehat{n}}.$$

(35)

Since the surface normal vector dependence can be expressed as $\widehat{n}=\widehat{n}(x,y)$ and by the application of the chain rule, we can derive:

$$\dot{\widehat{n}}={\widehat{n}}_x\dot{x}+{\widehat{n}}_y\dot{y},$$

(36)

and furthermore

$$\omega =\left[\begin{array}{cc}\widehat{n}\times {\widehat{n}}_x& \widehat{n}\times {\widehat{n}}_y\end{array}\right]\left[\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right],$$

(37)

where $\dot{x}$ and $\dot{y}$ are linear velocities in the $x$- and $y$- Cartesian coordinates, respectively. Now we can recall the Weingarten equations (16)-(19) and after some manipulation one can derive:

$$\omega =\frac{\left[\begin{array}{cc}{r}_x& r_y\end{array}\right]}{\sqrt{EG-F^2}}\left[\begin{array}{cc}M& N\\ -L& -M\end{array}\right]\left[\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right]=\frac{\left[\begin{array}{cc}{r}_x& r_y\end{array}\right]}{\sqrt{EG-F^2}}(A\cdot {\rm I}\text{}{\rm I})\left[\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right],$$

(38)

where ${\rm I}{\rm I}$ is the second fundamental form matrix (15), and $A$ is defined as:

$$A=\left[\begin{array}{cc}0& 1\\ -1& 0\end{array}\right].$$

(39)

The angular velocity vector $\omega$ is expressed by (38) as a linear combination of the tangent vectors $r_x$ and $r_y$. Therefore, it belongs to the tangent plane of the workpiece surface, and it will not produce tool rotation about itself. Furthermore, the tool angular velocity vector $\omega$ will be synchronized with the linear velocity in the xy-plane ${\left[\begin{array}{cc}\dot{x}& \dot{y}\end{array}\right]}^T$ by the matrix $S_{\omega }$:

$$S_{\omega }=\frac{\left[\begin{array}{cc}{r}_x& r_y\end{array}\right]}{\sqrt{EG-F^2}}\Pi =\frac{S_v\Pi }{\sqrt{EG-F^2}},$$

(40)

where we denote $\Pi =A\cdot {\rm I}\text{}{\rm I}$. Thus, we can express:

$$\omega =S_{\omega }\left[\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right],$$

(41)

where $S_{\omega }\in {ℝ}^{3\times 2}$ and its rank is $\mathrm{rank}(S_{\omega })\le 2$.

In the following, we rewrite the velocity constraints described by equations (33) and (41) such that they read as:

$$v_C=S_v\left[\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right],$$

(42)

and

$${\omega }_C=S_{\omega }\left[\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right],$$

(43)

where $v_C$ and ${\omega }_C$ represent linear velocity and angular velocity constrained to the surface, respectively. Note that $v_C$ (${\omega }_C$) belongs to the range-space of matrix $S_v$ ($S_{\omega }$). Furthermore, from (42) we can express

$$\left[\begin{array}{c}{\dot{x}}_C\\ {\dot{y}}_C\end{array}\right]=S_v^{†}{v}_C,$$

(44)

where we consider $v_C={\left[{\dot{x}}_C,{\dot{y}}_C,{\dot{z}}_C\right]}^T$, and the operator ${(.)}^{†}$ stands for the Moore-Penrose pseudoinverse of a matrix, such that $S_v^{†}={(S_v^T{S}_v^{ })}^{-1}{S}_v^T$ is the left pseudo-inverse of $S_v$. If we insert the last result into (43). It yields:

$${\omega }_C=S_{\omega }{S}_v^{†}{v}_C=S_C{v}_C.$$

(45)

Here, $S_C$ represents the so-called surface velocity constrained matrix, that can be expressed as:

$$S_C=S_{\omega }{S}_v^{†}=\frac{S_v\Pi S_v^{†}}{\sqrt{EG-F^2}}.$$

(46)

Note that $S_C\in {ℝ}^{3\times 3}$ and its rank is $\mathrm{rank}(S_C)\le 2$. It maps the linear velocity vector from the surface tangent plane to the angular velocity vector on the surface tangent plane.

3.3. The constrained robot kinematics model

The velocity kinematics play a central role in the robot manipulability analysis. In this Section, we formulate the problem as seeking the inverse velocity kinematics solution of (22)-(23) considering the velocity constraints given by synchronization expressions for linear and angular velocity (42)-(46), which constrain the motion of the end-effector on the path embedded on the workpiece surface, such that the end-effector maintains normal orientation. If, furthermore, we consider the joint velocity limitation given by (24) , then it will yield the constrained manipulability ellipsoid of a reduced dimension, as we will show later.

Similarly, as in [28], the Jacobian matrix is partitioned into two submatrices, i.e., the translational submatrix $J_T\in {ℝ}^{3\times n}$ and the rotational submatrix $J_R\in {ℝ}^{3\times n}$:

$$J=\left[\begin{array}{c}{J}_T\\ J_R\end{array}\right],$$

(47)

such that $J_T\dot{\theta }=v$ and $J_R\dot{\theta }=\omega$, where $v$ and $\omega$ are the linear velocity vector and the angular velocity vector in the robot operational space, respectively. In the following, we derive a solution of the velocity inverse kinematics considering the partitioned matrix(47). Baksalaray [56,57] discussed the problem of the Moore-Penrose inverse of a partitioned matrix. The application of such a solution formula in the case of the Jacobian (47) can be written by the augmented Jacobian submatrices that incorporate additional restrictions [28,41,58]:

$${\tilde{J}}_T=J_T{P}_R^n.$$

(48)

$${\tilde{J}}_R=J_R{P}_T^n.$$

(49)

where $P_T^n=\left(I-J_T^{†}{J}_T\right)$ and $P_R^n=\left(I-J_R^{†}{J}_R\right)$ are orthogonal projection matrices into the translational and rotational null space, respectively. Here, $I$ is the identity matrix of a proper dimension. Note that the definition of ${\tilde{J}}_T$ (${\tilde{J}}_R$) imposes such a constraint that the angular (linear) velocity is additionally zero, and thus it can be called a strong translational (rotational) Jacobian matrix [41]. The inverse of the partitioned squared Jacobian (n = 6, $J\in {ℝ}^{6\times 6}$) can now be expressed as the block matrix:

$$J^{-1}={\left[\begin{array}{c}{J}_T\\ J_R\end{array}\right]}^{-1}=\left[\begin{array}{cc}{\tilde{J}}_T^{†}& {\tilde{J}}_R^{†}\end{array}\right],$$

(50)

where ${\tilde{J}}_T^{†}\in {ℝ}^{n\times 3}$ and ${\tilde{J}}_R^{†}\in {ℝ}^{n\times 3}$. We can apply the inverse Jacobian in the block form (50) in the computation of the inverse velocity kinematics solution as

$$\dot{\theta }=J^{-1}\left[\begin{array}{c}v\\ \omega \end{array}\right]=\left[\begin{array}{cc}{\tilde{J}}_T^{†}& {\tilde{J}}_R^{†}\end{array}\right]\left[\begin{array}{c}v\\ \omega \end{array}\right]={\tilde{J}}_T^{†}v+{\tilde{J}}_R^{†}\omega .$$

(51)

The solution consists of two joint velocity components: the first component ${\tilde{J}}_T^{†}v$ is responsible for the end-effector translational motion without affecting its rotation, since it belongs to the null-space of $J_R$, whereas the second component ${\tilde{J}}_R^{†}v$ is responsible solely for the end-effector rotational motion with zero contribution to its translation, since it belongs to the null-space of $J_T$. Indeed, if we apply the solution (51) in the forward velocity kinematics mapping (22), it yields the desired result:

$$\begin{array}{c}J\dot{\theta }=\left[\begin{array}{c}{J}_T\\ J_R\end{array}\right]({\tilde{J}}_T^{†}v+{\tilde{J}}_R^{†}\omega )=\left[\begin{array}{c}{J}_T({(J_T{P}_R)}^{†}v+{(J_R{P}_T)}^{†}\omega )\\ J_R({(J_T{P}_R)}^{†}v+{(J_R{P}_T)}^{†}\omega )\end{array}\right]=\left[\begin{array}{c}{J}_T(P_R{(J_T{P}_R)}^{†}v+P_T{(J_R{P}_T)}^{†}\omega )\\ J_R(P_R{(J_T{P}_R)}^{†}v+P_T{(J_R{P}_T)}^{†}\omega )\end{array}\right]=\\ =\left[\begin{array}{c}{J}_T{P}_R{(J_T{P}_R)}^{†}v+J_T{P}_T{(J_R{P}_T)}^{†}\omega \\ J_R{P}_R{(J_T{P}_R)}^{†}v+J_R{P}_T{(J_R{P}_T)}^{†}\omega \end{array}\right]=\left[\begin{array}{c}Iv+0\omega \\ 0v+I\omega \end{array}\right]=\left[\begin{array}{c}v\\ \omega \end{array}\right]\end{array}.$$

(52)

Here, we exploit the well-known property of the Moore-Penrose pseudo inverse with the orthogonal projector, i.e. ${(J_{T/R}{P}_{R/T}^n)}^{†}=P_{R/T}^n{(J_{T/R}{P}_{R/T}^n)}^{†}$ [54,59]; consequently, $J_{T/R}{\tilde{J}}_{T/R}^{†}=I$ and $J_{T/R}{\tilde{J}}_{R/T}^{†}=0$.

In the following, we replace $v$ and $\omega$ in (51) with $v_C$ and ${\omega }_C$, respectively, and consider the expressions (42) and (43). It yields:

$$\dot{\theta }={\tilde{J}}_T^{†}{v}_C+{\tilde{J}}_R^{†}{\omega }_C=({\tilde{J}}_T^{†}{S}_v+{\tilde{J}}_R^{†}{S}_{\omega })\left[\begin{array}{c}\dot{x}\\ \dot{y}\end{array}\right].$$

(53)

Furthermore, we consider (44) , and then (53) it can be rewritten as:

$${\dot{\theta }}_C=({\tilde{J}}_T^{†}{S}_v{S}_v^{†}+{\tilde{J}}_R^{†}{S}_{\omega }{S}_v^{†})v_C=J_C^{\#}{v}_C,$$

(54)

where ${\dot{\theta }}_C$ stands for such joint velocities, which provide robot tool motion constrained on the workpiece surface, and $J_C^{\#}\in {ℝ}^{n\times 3}$ (with $n=6$) denotes the so-called inverse Jacobian matrix of the constrained kinematics,

$$J_C^{\#}={\tilde{J}}_T^{†}{S}_v{S}_v^{†}+{\tilde{J}}_R^{†}{S}_{\omega }{S}_v^{†}.$$

(55)

Note that $S_v{S}_v^{+}$ represents the orthogonal projector to the range space of matrix $S_v$ (with the following properties ${(S_v{S}_v^{+})}^2={(S_v{S}_v^{+})}^T=S_v{S}_v^{+}$), such that $v_C=S_v{S}_v^{+}v$, and from (45) it follows directly that $S_{\omega }{S}_v^{+}$ maps any $v_C$ to the range space of matrix $S_{\omega }$. Furthermore, if we consider the orthogonal decomposition of the operational space velocity vector $v$, such that:

$$v=v_C+v_C^{\perp },$$

(56)

where $v_C^{\perp }=(I-S_v{S}_v^{+})v$ is orthogonal to $v_C$, i.e. $v_C^{\perp }\perp v_C$, then it is easy to show that (54) it can also be written as:

$${\dot{\theta }}_C=({\tilde{J}}_T^{†}{S}_v{S}_v^{†}+{\tilde{J}}_R^{†}{S}_{\omega }{S}_v^{†})v=J_C^{\#}v$$

(57)

since matrix $J_C^{\#}$ maps the velocity component vector $v_C^{\perp }$ to zero:

$$\begin{array}{c}{J}_C^{\#}{v}_C^{\perp }=({\tilde{J}}_T^{†}{S}_v{S}_v^{†}+{\tilde{J}}_R^{†}{S}_{\omega }{S}_v^{†})(I-S_v{S}_v^{†})v=({\tilde{J}}_T^{†}{S}_v{S}_v^{†}(I-S_v{S}_v^{†})+{\tilde{J}}_R^{†}{S}_{\omega }{S}_v^{†}(I-S_v{S}_v^{†}))v=\\ =({\tilde{J}}_T^{†}(S_v{S}_v^{†}-(S_v{S}_v^{†})2)+{\tilde{J}}_R^{†}{S}_{\omega }(S_v^{†}-S_v^{†}{S}_v{S}_v^{†}))v=({\tilde{J}}_T^{†}(S_v{S}_v^{†}-S_v{S}_v^{†})+{\tilde{J}}_R^{†}{S}_{\omega }(S_v^{†}-S_v^{†}))v=0\end{array}.$$

(58)

Finally, we show that the derived solution above provides such velocities in the operational robot space, which are constrained to the workpiece tangent plane. It is rather straightforward to derive that the forward velocity kinematics mapping (22) on ${\dot{\theta }}_C$ results as

$$J{\dot{\theta }}_C=J(J_C^{\#}{v}_C)=\left[\begin{array}{c}{J}_T\\ J_R\end{array}\right]({\tilde{J}}_T^{†}{S}_v{S}_v^{†}+{\tilde{J}}_R^{†}{S}_{\omega }{S}_v^{†})v_C=\left[\begin{array}{c}{J}_T{\tilde{J}}_T^{†}{S}_v{S}_v^{†}+J_T{\tilde{J}}_R^{†}{S}_{\omega }{S}_v^{†}\\ J_R{\tilde{J}}_T^{†}{S}_v{S}_v^{†}+J_R{\tilde{J}}_R^{†}{S}_{\omega }{S}_v^{†}\end{array}\right]v_C=\left[\begin{array}{c}{S}_v{S}_v^{†}\\ S_{\omega }{S}_v^{†}\end{array}\right]v_C=\left[\begin{array}{c}{v}_C\\ {\omega }_C\end{array}\right]$$

(59)

which confirms the solution described by (54), (55) and (57) provides the desired constrained motion (42)-(43) in the robot operational space. Furthermore, if we consider that the inverse Jacobian of the constrained kinematics can also be represented as $J_C^{\#}=(J^{-1}S)S_v^{†}$, where $S^T=\left[\begin{array}{cc}{S}_v^T& S_{\omega }^T\end{array}\right]$, it can also be shown easily that its rank equals $\mathrm{rank}(J_C^{\#})=2$; on the basis of the Sylvester rank inequality formula [60] it is possible to derive$\mathrm{rank}(J^{-1}S)=2$ and in the following also $\mathrm{rank}((J^{-1}S)S_v^{†})=2$, since the matrix $S_v^{†}$ is of full rank, i.e. $\mathrm{rank}(S_v^{†})=2$.

3.4. The manipulability ellipsoid of the constrained motion

Now, we can construct the manipulability ellipsoid of the constrained motion following a procedure that is similar to the steps determined by (24)-(29). Firstly, we consider the constraint in the joint space given by the condition (24) , that is derived as:

$${\dot{\theta }}^T\dot{\theta }=v^T(J_C^{\#}{ }^T{J}_C^{\#})v\le 1,$$

(60)

where $J_C^{\#}$ is defined by (55). In the next step, we can factorize matrix $J_C^{\#}$ by SVD into

$$J_C^{\#}=U\Sigma V^T.$$

(61)

However, due to the fact that the matrix is tall and its rank is deficient $r=\mathrm{rank}(J_C^{\#})=2$, the factorized matrices $U$, $\Sigma$ and $V$ have the following structure:

$$U=\left[\begin{array}{cc}{U}_1& U_2\end{array}\right],\Sigma =\left[\begin{array}{cc}{\Sigma }_1& 0\\ 0& 0\end{array}\right],\mathrm{and}V=\left[\begin{array}{cc}{V}_1& V_2\end{array}\right]$$

(62)

respectively. Here, the columns of $U_1$ represent an orthonormal base of the column space $Col(J_C^{\#})$, whereas the columns of $U_2$ represent an orthonormal base of the left null space $Null(J_C^{\#}{ }^T)$. Furthermore, the columns of $V_1$ represent an orthonormal base of the row space $Row(J_C^{\#})$, whereas the columns of $V_2$ represent an orthonormal base of the null space $Null(J_C^{\#})$. The (sub)matrices` dimensions are $\dim (U_1)=n\times r$, $\dim (U_2)=n\times (n-r)$, $\dim (V_1)=3\times r$, $\dim (V_2)=3\times (3-r)$ and $\dim ({\Sigma }_1)=r\times r$, respectively. Note that $U$and $V$ are unitary matrices with orthonormal columns, and ${\Sigma }_1$ is a diagonal square matrix with the positive singular values ${\sigma }_1\ge \mathrm{...}\ge {\sigma }_r>0$ on the diagonal. We can derive:

$$J_C^{\#}=\left[\begin{array}{cc}{U}_1& U_2\end{array}\right]\Sigma =\left[\begin{array}{cc}{\Sigma }_1& 0\\ 0& 0\end{array}\right]\left[\begin{array}{c}{V}_1^T\\ V_2^T\end{array}\right]=U_1{\Sigma }_1{V}_1^T,$$

(63)

and further

$$J_C^{\#}{ }^T{J}_C^{\#}={(U_1{\Sigma }_1{V}_1^T)}^T(U_1{\Sigma }_1{V}_1^T)=(V_1{\Sigma }_1{U}_1^T)(U_1{\Sigma }_1{V}_1^T)=V_1{\Sigma }_1{U}_1^T{U}_1{\Sigma }_1{V}_1^T=V_1{\Sigma }_1^2{V}_1^T.$$

(64)

Thus, (60) it can be rewritten as

$$v^T{V}_1{\Sigma }_1^2{V}_1^{T}v\le 1.$$

(65)

If we introduce a new variable of reduced order $w=V_1^{T}v$, then it yields the ellipsoid,

$$w_{ }^T{\Sigma }_1^{2}w\le 1$$

(66)

that can be rewritten as

$$\frac{w_1^2}{{\sigma }_1^{-2}}+\frac{w_2^2}{{\sigma }_2^{-2}}\le 1,$$

(67)

since $r=2$, and thus $w_{ }^T=\left[\begin{array}{cc}{w}_1^{ }& w_2^{ }\end{array}\right]$, $\mathrm{diag}({\Sigma }_1)=\left[\begin{array}{cc}{\sigma }_1& {\sigma }_2\end{array}\right]$. Obviously, the constrained manipulability ellipsoid is reduced to an ellipse in the tangent plane of the workpiece surface.

4. Numerical Experiments and Discussion

In this Section we show the results of some numerical experiments that we performed in order to verify the proposed approach of the robot manipulability evaluation in the case of the robot tool motion constrained by the workpiece surface. In the numerical experiments we applied a collaborative robot arm and two workpiece models with different surface characteristics. In the following, we firstly introduce the experimental setup. Then, we show the ellipsoids related to the velocity synchronization on the workpiece surface, which can be derived from (45)-(46). Next, we present the manipulability ellipsoids of the robot motion with the tool constrained to the surface geometry derived from (54)-(55), and it is described further in Section 3.4. Finally, we verify the drawn manipulability ellipsoids by the calculation of the robot tool velocity in selected directions. The numerical experiments were performed by the MATLAB computing environment [61].

4.1. Experimental setup

The robotic setup is shown in Figure 1. We applied a collaborative robot arm UR5 with a reach radius of 850 (mm); the kinematic parameters of the robot arm are available in reference [62]. The reference robot frame is shown by the red-green-blue axes located at the robot base. The MHP tool length was modeled as 284.5 (mm). The workpiece size was 250 x 250 x 75 (mm), and was placed with its center of the bottom face at the point [0, -525, -97] (mm) of the reference robot frame. In our numerical experiments, we applied two workpiece surface models, i.e., Workpiece surface model No.1 (WP1) and Workpiece surface model No.2 (WP2), which are shown by Figure 2 and Figure 3, respectively. WP1 had a relatively monotone curvature, while WP2 had a more vivid curvature. The colormap shown on the workpiece surface is coded by its height.

4.2. Surface velocity ellipsoids

The ellipsoids related to the velocity synchronization on the workpiece surface are presented by Figure 4, Figure 5 and Figure 6. They are developed in the selected points on the workpieces’ surface, and represent the omnidirectional angular velocity (${\omega }_C$) performance of a particle moving along the path constrained on the workpiece surface under the condition $‖v_C‖\le 1$. The surface velocity synchronization ellipsoids, which can be derived from (45)-(46), are shown firstly by Figure 4 for WP1 and Figure 5 for WP2. The ellipsoids are clearly of reduced dimension, such that they convert to ellipses and appear in the tangent planes attached to the selected points on the surface. The shown ellipses were all scaled by the factor of 0.00075. The angular velocity ellipses were of different sizes and eccentricity. Thus, the angular velocity performance depended both on the selected direction of the linear movement, and on the geometry of the close neighborhood of the selected point on the workpiece surface. If a big ellipse similar to a circle was drawn, then a relatively high angular velocity would occur in all directions of the linear movement. In the case of an eccentric ellipse, the highest angular velocity would occur in the direction of the major ellipse axis, and the lowest in the direction of the minor ellipse axis. The angular velocity performance depends highly on the surface slope (the surface gradient of the first order), its gradient (the surface gradient of the second order), and, thus, on the surface curvature.

The same angular velocity ellipses are again shown in Figure 6. However, in this case, the surface colormap is coded by its Mean curvature measure. A red/blue color means a high/low curvature measure value. It can be observed clearly that the drawn ellipses’ characteristics correspond highly to the surface curvature measure at the selected points on both WP1 and WP2: at the surface points with low value of the curvature measure, we can notice small ellipses, whereas at the surface points with high value of the surface curvature measure, we can see big ellipses; in the case of the medium value curvature points, we can observe medium size ellipses, or ellipses with a longer major axis (and a medium average axes` length), that again match with the average curvature at the surface point. We can conclude that the angular velocity ellipses shown obviously correlate with the surface curvature characteristics.

4.3. Robot manipulability ellipsoids of the constrained kinematics

The ellipsoids related to the robot translational velocity along the path constrained on the workpiece with a curved surface and with the robot tool orientation normal to the surface are presented by Figure 7, Figure 8, Figure 9 and Figure 10. They express the robot linear velocity performance at the positions in the robot operational space determined by the selected points on the workpiece surface under the condition of unity Euclidian norm of the robot joint velocities $‖\dot{q}‖\le 1$. The constrained manipulability ellipsoids, which can be derived from (54)-(55), are shown firstly by Figure 7 for WP1 and Figure 8 for WP2. The ellipsoids are again clearly of reduced dimensions, such that they convert to ellipses and appear in the tangent planes attached to the selected points on the surface. The shown ellipses were all scaled by a factor of 0.05.

At the bottom of WP1, all around in Figure 7 (b), we can see ellipses with different sizes and orientations. Most of them are of medium size; some ellipses at the bottom side of the picture are longer and thinner, that determines a larger linear velocity performance in the direction of the ellipses’ major axes, and poorer in the direction of the minor axes; however, on the upper side of the picture the ellipses are of small size, which is connected to poor linear velocity performance in all directions. Along the slope of WP1 we can see many rather long ellipses, with a major axis directed towards the top of the workpiece; it means that a straightway motion from the bottom to the top of the surface would be allowed at relatively high speed. Finally, on the top of the workpiece we can observe very small velocity ellipses, since the surface in this area changed rapidly.

The velocity manipulability ellipsoids on WP2, with more vivid curvature surface in Figure 8, are with different sizes, orientations, and widths. Their shape depends on the surface characteristics at the specific point and the specific robot configuration required to provide the robot tool position and orientation that comply with the imposed surface constraints. Some big/medium size and only slightly eccentric ellipses, which can be observed mostly in the bottom half of the surface presented by Figure 8 (b), allow for high/medium linear omnidirectional velocity performance along the tangent planes attached to certain surface points. Thin ellipses on highs at the right side and in lows at the left side provide better velocity performance in certain directions determined by the ellipses’ major axis. Many small size ellipses at different points all around the workpiece surface show only poor linear velocity performance.

A more detailed insight into the linear velocity performance is presented by Figure 9, where the surface colormap is coded by the mean (average) value of the ellipses’ minor and major axes` lengths at all surface definition points. The red/blue color means roughly good/bad velocity manipulability. Figures 9 (a)/(b) show the surface constrained robot velocity manipulability for WP1/WP2. These characteristics can be used as a robot tool velocity regulation map on the surface at a certain workpiece location, since they incorporate not only robot kinematics performance, but also the workpiece surface curvature feature.

Finally, we verified the presented velocity ellipses for both WP1 and WP2. In the verification procedure we simply calculated the robot tool linear velocity of maximal possible speed that can be achieved under condition $‖\dot{q}‖\le 1$ during the motion along the path from the selected point to a certain neighboring point. In every ellipse we selected 8 neighboring points around the ellipse center point. In order to determine the maximum linear speed, we applied the following formula:

$$V_{\max }=\frac{1}{{‖{\tilde{J}}_T^{+}{u}_T+h^{-1}{\tilde{J}}_R^{+}{u}_R‖}_2},$$

(68)

where $u_T$ and $u_R$ are linear and angular directions along the path, whereas the parameter $h$ is related to the surface shape/curvature in our case, and can be calculated as $h=s/\phi$ where $s$ and $\phi$ stand for linear path length and rotational path length (angle), respectively. The formula above is similar to the DTF method formula [28], however, please note that the DTF method employs the infinity norm in the denominator, whereas we used the Euclidean L2 norm. The verification results for WP1/WP2 are shown by Figures 10 (a)/(b). Here, the surface colormap is coded by the surface height. Again, we can see the before mentioned manipulability ellipses, as in Figure 7 and Figure 8. Beside these, there are the verification linear velocity vectors drawn in the bright pale pink color. They were all scaled by the factor of 0.05, which is of the same value as the manipulability ellipses` scaling factor. We can observe that the vector lines drawn from the ellipses’ center point end at the ellipses’ edge. Thus, they correspond to the velocities, which could be determined from the manipulability ellipses in the verification directions, and therefore confirm the correctness of the presented manipulability ellipses of the robot kinematics imposed by the workpiece surface constraints.

5. Conclusions

In this paper, we discussed the robot kinematics` capability and velocity performance in the case of workpiece machining. We consider such machining tasks as the hammer-peening workpiece surface polishing, in which the robot tool is required to follow a prescribed path on the workpiece with zero exerted force in the path direction that enforces certain kinematics constraints to the robot`s motion. In the presented study, we have considered a workpiece with complex surface characteristics, to derive the robot`s constrained kinematics description. Thus, we introduced a so-called surface velocity constrained matrix, which describes the synchronization of the linear velocity and the angular velocity constrained on the surface tangent plane. The surface velocity synchronization was then applied in the derivation of the constrained robot kinematics. The derived Jacobian kinematics matrix with incorporated workpiece surface constraints was, in the following, used for the ellipsoid-based manipulability analysis. We have shown that, in this case, the general non-homogeneous 6-dimensional manipulability ellipsoids, which mix different units of linear velocity and angular velocity, can be reduced to the 2-dimensional ellipses in the homogeneous space of translational motion. The results from the theoretical development have been proven by numerical experiments. In the future, the presented approach of the surface-based constrained manipulability can be employed in the kinematic optimization of the robot-workpiece position and planning of the optimal robot machining path. Furthermore, future experimental validation may include real world manufacturing examples in a real robotic cell with a collaborative robot.