Kinematic of 4SPRR-SPR parallel robots based on VQTAM

Abstract

4SPRR-SPR parallel robot is a novel closed-loop mechanism. The research of kinematics is the basis of real-time and robust robot control. This paper aims to proposing a method to address a surrogate model of forward kinematics for PMs (Parallel Mechanisms). Herein, the forward kinematics model is derived by training the VQTAM (Vector-Quantified Temporal Associative Memory) network, which originates from a SOM (Self-Organized Mapping). During the processes of training, testing and estimating this neural network, the priority K-means tree search algorithm is utilized, thus improving the training efficacy. Furthermore, LLR (Local Linear Regression), LWR (Local Weighted Linear Regression) and LLE (Local Linear Embedding) algorithms are respectively combined with VQTAM to get three improvement algorithms, aiming to further optimize the prediction accuracy of the networks. To speed up solving the least squared equation in the three algorithms, SVD (Singular Value Decomposition) is introduced. Finally, Data from inverse kinematics by geometric method is obtained, which is for constructing and validating the VQTAM neural network. Results show that the prediction effect of LLE algorithm is better than others, which could be a potential surrogate model to estimate the output of forward kinematics.

Introduction

Parallel robots offer noticeable advantages in machining applications, especially in the field of heavy equipment manufacture, due to their high rigid-weight ratio and load capacity. Recently, the machine tool industry has discovered the potential advantages of PMs and many parallel machine tools have been developed based on either the 6-DOF (degrees-of-freedom) parallel mechanisms or 5-DOF PMs. ¹

Parallel robots or PMs with 5 or 6 DOF are often used as motion simulators and crane devices. Stewart introduced a 6-DOF PM, which is now popularly known as Stewart platform. ² It is commonly used in flight or driving simulators, vibration isolation platform, space docking mechanisms, and space telescope.^3,4 Hereafter the lower-mobility (2-5 DOFs) PMs have emerged and have become popular owing to their simple structure, low kinematic coupling, low cost and convenient control.^5–7 Xie et al., ⁵ Huang et al. ⁶ and Joshi et al. ⁷ studied the type synthesis for different types of lower-mobility PM. Therefore the 2UPR/SPR PM ⁸ and 3UPS/UP PM^9,10 have been used as position adjusters and space telescope (R, P, U and S stand for revolute, prismatic, universal and spherical joints respectively). As far as 5-DOF PMs are concerned, researchers have mainly worked on the type synthesis.^6,11–15 Through the type synthesis, a variety of PMs are obtained.^16,17 For instance, M. Garcia-Murillo et al. ¹⁸ proposed a 2(3-RRPS) parallel robot and analyzed its workspace. Nicholas Baron, et al. ¹⁹ displayed a family of kinematically redundant planar parallel robots, and presented a geometric method of singularity avoidance. In addition, the cable-driven manipulators are a class of PMs, where the end-effector is displaced by means of servo-actuated extendable cables. ²⁰ They have high ratio of payload to robot weight and a much larger workspace which is appropriate for camera positioning in sports stadiums and many others. ²¹ Arian Bahrami et al. studied the optimal design of various forms of the cable-driven manipulators.^22,23 In summary, the structures of PMs are various and have been applied to many fields. However, there is still no executable method to solve the kinematics, especially the forward kinematics.

The 5SPRR PM is a 6-DOF PM, which has been used as a machining tool.^24–26 The structural form of 4SPRR-SPR parallel robot, which has great potential to transform into a new type of mobile machining robot, is derived from the 5SPRR PM. To the best knowledge of the author, up to now, Research, such as its kinematic, workspace analysis et al. has not been conducted yet.

Parallel manipulators do not have formulaic kinematic solutions. Therefore, several studies are focused on methods to simplify the solution of robot kinematics. Intelligent algorithm or machine learning is introduced into kinematics solution of robot. Simpler mapping models of kinematic input and output are established to replace the complex analytical solution. Evolutionary computational methods such as genetic algorithm (GA),^27–31 differential evolution (DE), and particle swarm optimization (PSO) ³² are used to optimize the solutions by formulating the problem in the form of an objective or a fitness function. Furthermore, approximation methods based on artificial neural networks and polynomial regression models^33–37 can also be used to estimate forward kinematic. Robot kinematics can be regarded as a non-linear autoregressive model with exogenous input (NARX). In reference, ³⁸ the forward kinematics of parallel manipulators is solved by combining the wavelet network and the NARX method. However, the validity of this method is still dependent on parameter initialization and back propagation learning process. The simplification of complex kinematics model often leads to decreasing the accuracy. To overcome this problem, efficient and accurate kinematics identification is needed. ³⁹

The VQTAM algorithm does not depend on the adjustment of super parameters, but on SOM to achieve topology preservation. It can be applied to the mapping from input space to output space, as an extension of SOM. The mathematical foundation of VQTAM is non-linear manifold embedding. This method is suitable for time series modeling and prediction of non-linear systems. Depending on the self-organizing mapping structure, the discrete mapping relationship can effectively realize output estimation according to exogenous input.^39–41 Since the input space and output space of VQTAM are discrete, a network with more neurons is required to achieve high-precision prediction. In order to improve the computational efficiency, it is necessary to reduce the number of VQTAM neurons, which decreases the efficiency. Improved local linear algorithms of VQTAM are proposed to solve the inverse kinematics of a 6-DOF robot efficiently and accurately.

In this paper, a novel parallel robot that could be used as mobile machining robot is introduced. The architecture and general kinematic properties of the 4SPRR-SPR parallel robots are outlined, and the study on its kinematic follows.

Research methodology

4SPRR-SPR parallel robot is a novel closed-loop mechanism. Its inverse kinematics can be achieved using a geometric method. The rod length can be obtained quickly through the known pose of end effector. But the forward kinematics of parallel robot is very complex. It’s a non-linear system. Generally, it is difficult to accurately determine the forward kinematics.

Thus, the length of drive limbs corresponding to different poses can be obtained using inverse kinematics conveniently. The data obtained from the inverse kinematics can be used to train the surrogate model of the forward kinematics. In this paper, VQTAM neural network is used to replace the complex calculation of forward kinematics. The research content and structure of this paper are shown in Figure 1.

Figure 1.

Research content and structure.

In Section 2, the inverse kinematics is studied. Using the inverse kinematics algorithm, the length of drive limbs corresponding to poses is calculated as the sample set.

In Section 3, the non-linear system identification method based on VQTAM is studied, and the training method of VQTAM neural network is established. In the process of VQTAM training and testing, it is core to search and activate neuron nodes. Therefore, the priority search K-means tree algorithm is introduced for this step.

In Section 4, the improved method of local linearization of input and output mapping based on VQTAM neural network is studied.

Kinematic of robotic manipulator

4SPRR-SPR parallel robots

The 4SPRR-SPR parallel robot is one of the new 5-DOF PMs, and its structure is shown in Figure 2. The four limbs of the PM consist of one spherical joint, one moving pair and two rotating pairs, and the other one lack for a rotating pair at one end of the limb. The rotating pairs connected with an end effector are coaxial. The spherical joints connected with the static platform are not arranged in the same plane but are distributed in different positions in space. This provides more possibilities for the optimization of the parallel mechanism.

Figure 2.

The structure of 4SPRR-SPR parallel robot.

The arrangements of the kinematic pairs of the 4SPRR-SPR parallel robot are similar, except that the end of the first limb is fixed to the moving platform. The kinematic pair arrangement of the i^th chain is shown in Figure 3 (i = 2,3,…5).

Figure 3.

ith limb structure diagram.

A_i is the center point of the universal joint connected by the ith limb and the static platform, B_i and C_i are the two revolving pair centers of the limb, and D is the reference point of the position of the end effector. The screw of the center line S of each pair is used to represent the motion of the pair. ¹² S _i ¹ and S _i ² are revolting joints; S _i ³ is a prismatic joint, and S _i ⁴, S _i ⁵, and S _i ⁶ constitute the rotation of spherical joint in three directions.

For the ith limb, the coordinate system O_i-x_iy_iz_i at E_i is established. E_i is the vertical foot from A_i to DC_i. The screws represented in coordinated O_i-x_iy_iz_i are

$$\{\begin{array}{c}{\mathit{S}}_i^1=(0,0,1;0,0,0)\hfill \\ \multicolumn{1}{c}{{\mathit{S}}_i^2=(0,1,0;m_i,0,c_i)}\\ \multicolumn{1}{c}{{\mathit{S}}_i^3=(0,0,0;n_i,0,m_i)}\\ \multicolumn{1}{c}{{\mathit{S}}_i^4=(0,1,0;0,0,n_i)}\\ \multicolumn{1}{c}{{\mathit{S}}_i^5=(1,0,0;0,0,0)}\\ \multicolumn{1}{c}{{\mathit{S}}_i^6=(0,0,1;0,-n_i,0)}\end{array}$$

(1)

Where b_i, c_i, m_i, and n_i represent the length of segment A_iB_i, B_iC_i, C_iE_i, and E_iA_i respectively.

The six screws of the i^th limb are linearly independent (i = 2,3,…5). The screws of the i^th limb constitute a 6-system of screws, which provide no constraint to the end effector.

The screw system of the first limb differs from those of others, because the end of the first limb is fixed to the moving platform. The screws represented in coordinated O₁-x₁y₁z₁ are

$$\{\begin{array}{c}{\mathit{S}}_1^2=(0,1,0;m_1,0,c_1)\hfill \\ \multicolumn{1}{c}{{\mathit{S}}_1^3=(0,0,0;n_1,0,m_1)}\\ \multicolumn{1}{c}{{\mathit{S}}_1^4=(0,1,0;0,0,n_1)}\\ \multicolumn{1}{c}{{\mathit{S}}_1^5=(1,0,0;0,0,0)}\\ \multicolumn{1}{c}{{\mathit{S}}_1^6=(0,0,1;0,-n_1,0)}\end{array}$$

(2)

The reciprocal screw of the 5-system screws can be written as

$${\mathit{S}}_1^r=(0,1,0;0,0,n_1)$$

(3)

Therefore, five limbs provide one constraint to the end effector. It can be confirmed that 4SPRR-SPR parallel robots have 5 DOFs. The motion screw system which reflects the DOF property of the end-effector can be obtained by finding the secondary inverse screw system of the end-effector.

$$\begin{array}{c}{\mathit{S}}_1^{\mathrm{pm}}={\mathit{S}}_1^{\mathrm{rr}}=(1,0,0;0,0,0)\hfill \\ {\mathit{S}}_2^{\mathrm{pm}}={\mathit{S}}_2^{\mathrm{rr}}=(0,1,0;0,0,0)\hfill \\ {\mathit{S}}_3^{\mathrm{pm}}={\mathit{S}}_3^{\mathrm{rr}}=(0,0,0;1,0,0)\hfill \\ {\mathit{S}}_4^{\mathrm{pm}}={\mathit{S}}_4^{\mathrm{rr}}=(0,0,0;0,0,1)\hfill \\ {\mathit{S}}_5^{\mathrm{pm}}={\mathit{S}}_5^{\mathrm{rr}}=(0,0,1;0,-n_1,0)\hfill \end{array}$$

(4)

Formula (4) shows that the end effector can rotate around the X and Y axes and around the line (0,0,1;0,-n₁,0), which passes through point A₁ and parallels to the Z axis. The end effector can also moving along the X and Z-axes.

Inverse kinematics of 4SPRR-SPR parallel robot

The prismatic joints of the five limbs constitute the drive of a 4SPRR-SPR parallel robot. Therefore, the inverse position solution of the parallel robot can be regarded as the calculation of the rods A_iB_i about the length l_i using the coordinates of the known points D and direction vectors S ¹ (i = 1,2,…,5). The directions of end effector are the same for the five limbs in unified coordinates.

The geometric condition of the parallel robot is determined using equations (1) and (2).

$$\{\begin{array}{c}{\mathit{c}}_i\perp {\mathit{S}}^1\hfill \\ \multicolumn{1}{c}{{\mathit{S}}_i^2\perp {\mathit{c}}_i}\\ \multicolumn{1}{c}{{\mathit{b}}_i\perp {\mathit{S}}_i^2}\\ \multicolumn{1}{c}{{\mathit{S}}^1\perp {\mathit{S}}_i^2}\end{array}$$

(5)

Where ${\mathit{b}}_i={\mathit{A}}_i{\mathit{B}}_i$ , ${\mathit{c}}_i={\mathit{B}}_i{\mathit{C}}_i$ , ${\mathit{d}}_i={\mathit{C}}_i{\mathit{D}}_i$

As shown in Figure 3, the i^th limb is considered for the analysis. E_i is the vertical foot from A_i to DC_i. From the geometric relation, the result is deduced using equations (6) to (8).

$${\mathit{E}}_i{\mathit{C}}_i=({\mathit{A}}_i{\mathit{C}}_i·{\mathit{S}}^1){\mathit{S}}^1$$

(6)

$${\mathit{A}}_i{\mathit{E}}_i={\mathit{A}}_i{\mathit{C}}_i-{\mathit{E}}_i{\mathit{C}}_i$$

(7)

$${\mathit{c}}_i=c_i\left(\frac{{\mathit{A}}_i{\mathit{E}}_i}{\left|{\mathit{A}}_i{\mathit{E}}_i\right|}\right)$$

(8)

The coordinates of points C_i and B_i can be obtained upon calculating c_i. Accordingly the input parameter l_i corresponding to the pose can be obtained.

The main structural parameters of 4SPRR-SPR parallel robots are the coordinates of points A_i (i = 1,2,…5) and rod lengths c_i and d_i. The above parameters can determine the corresponding parallel robot for inverse kinematics calculation.

The inverse kinematics solution of the robot is calculated considering the structural parameters shown in Table 1 as an example.

Table 1.

The structural parameters of the parallel robot.

Parameters	Values (mm)
A 1	(0, 2380.19, 1823.98)
A 2	(−790.54, 1011.96, 1877.15)
A 3	(790.54, 1011.96, 1877.15)
A 4	(−557.90, 2471.63, 2218.87)
A 5	(557.90, 2471.63, 2218.87)
d 1	440
d 2	520
d 3	600
d 4	1040
d 5	1120
ci	180

By determining the position and posture of the moving platform, the coordinate of point D and the direction vector S ¹ is found to be (0,1119.17,3641.27), and (0,0,-1). The length l_i of the rod A_iB_i, as seen in Table 2, is calculated by substituting these parameters in the aforementioned algorithm.

Table 2.

The calculations of parallel robot.

The rod length label	The rod length(mm)
l 1	1751.26
l 2	1386.05
l 3	1317.89
l 4	1338.79
l 5	1318.17

System identification of non-linear dynamical system using VQTAM

The forward kinematics of the parallel robot aims at solving the corresponding pose given the length of drive limbs. The pose of the end effector is regarded as output Q and the length of the drive limbs as input p . The forward kinematics problem of the robot can be transformed into the identification problem of a non-linear dynamic system. The time-discrete difference formula is established as shown in Formula (9).^39,41

Where Q (t) =[x(t),y(t),z(t); S ¹(t)], which is the pose of the end effector; p (t) = [l₁(t),l₂(t),l₃(t),l₄(t),l₅(t)], which is the length of drive limbs at t time. f(·) represents a non-linear function reflecting the characteristics of the system. The output Q at t+1 is determined by the previous n_q outputs and n_p inputs.

The purpose of forward kinematics is to recognize the mapping relation f(·) from input to output. For a non-linear problem without formulaic solutions, simplifying the model can guarantee efficiency and ensure favorable real-time performance in engineering applications.

$$\mathit{Q}(t+1)=f[\mathit{Q}(t),...,\mathit{Q}(t-n_q+1);\mathit{p}(t),...,\mathit{p}(t-n_u+1)]$$

(9)

$${\mathit{X}}^{\mathrm{in}}=[\mathit{Q}(t),...,\mathit{Q}(t-n_q+1);\mathit{p}(t),...,\mathit{p}(t-n_u+1)]$$

(10)

$${\mathit{X}}^{\mathrm{out}}=[\mathit{Q}(t+1)]$$

(11)

VQTAM is an extension of SOM, which can be used as an effective manner to realize non-linear mapping. The concept of non-linear manifold embedding can be used to establish the non-linear mapping. SOM is realized through competition and cooperation among feature neurons. Different neurons are activated for different inputs, and the activated neurons affect the output together with the neurons in their neighborhood. Therefore, time series data can be constructed as the training samples of VQTAM.

VQTAM uses the training samples to build a topological structure embedded in multi-dimensional data space, which can construct the mapping relationship between multi-dimensional input vectors and multi-dimensional output vectors. The VQTAM network consists of three layers: input space ω ⁱⁿ, output space ω ^out and lattice space. The establishment of the mapping process is shown in Figure 4. The definitions of X ⁱⁿ and X ^out are shown in equations (10) and (11).

Figure 4.

VQTAM spatial mapping.

The input space and output space are composed of weight vectors ω _i ⁱⁿ and ω _i ^out respectively. i is the index of neuron in lattice space, which reflects the location of neuron in lattice topology structure. The weight vectors ω _i ⁱⁿ and ω _i ^out have mapping relationships with neurons in lattice space. ω _i ⁱⁿ and ω _i ^out have the same dimensions as X ⁱⁿ and X ^out respectively. In the mapping process of VQTAM, the nearest weight vector ω _i* ⁱⁿ to X ⁱⁿ is found in the input space, which corresponds to the activated neuron in the lattice space. In this process, The Euclidean distance is used to measure the distance between two vectors. The definition of i* is shown in Formula (12), where A is the collection of all neuron indexes in lattice space. ⁴¹

$$i^{*}=\underset{i\in A}{\arg min}\{\Vert {\mathit{X}}^{\mathrm{in}}-{\mathit{\omega }}_i^{\mathrm{in}}\Vert \}$$

(12)

ω _i* ^out is obtained according to the neuron index i* and the output X^out is estimated, as shown in Formula (13).

$${\hat{\mathit{X}}}^{\mathrm{out}}={\mathit{\omega }}_{i*}^{\mathrm{out}}$$

(13)

Learning strategy of VQTAM

The main parameters of VQTAM are neuron weight vector ω _i ⁱⁿ, ω _i ^out and lattice topology. The initial lattice topology structure can be determined as a super parameter before training, so the learning strategy mainly aims at updating the weight vectors ω _i ⁱⁿ and ω _i ^out in the training process.

There is competition and cooperation among neurons in SOM. Different neuron and their neighborhoods are activated by input to participate in the learning process, and the weight vectors ω _i ⁱⁿ and ω _i ^out are updated. The index search of the active neurons in the lattice space is consistent with equation (12).

In order to improve the convergence rate in the learning process, the influence range parameters σ(t) of the activation neuron and the learning speed α(t) decrease exponentially with the learning epoch, as shown in equations (14) and (15).

$$\alpha (t)={\alpha }_0({\alpha }_T/{\alpha }_0{)}^{t/T}$$

(14)

$$\sigma (t)={\sigma }_0({\sigma }_T/{\sigma }_0{)}^{t/T}$$

(15)

Where t is the current learning epoch, and T is the total learning epoch. α₀ and σ₀ are the initial value, and α_T and σ_T are the parameter value of training epoch T.

The Gauss neighborhood function is used to determine the effect of input on the neighbors of activated neuron, as shown in equation (16).^39–41

$$h(i^{*},i;t)=\exp (\frac{-{\Vert i(t)-i^{*}(t)\Vert }^2}{{\sigma }^2(t)})$$

(16)

ω _i ⁱⁿ and ω _i ^out are updated, as shown in equations (17) and (18).

$${\mathit{\omega }}_i^{\mathrm{in}}\Leftarrow \alpha (t)h(i^{*},i;t)[{\mathit{X}}^{\mathrm{in}}-{\mathit{\omega }}_i^{\mathrm{in}}]$$

(17)

$${\mathit{\omega }}_i^{\mathrm{out}}\Leftarrow \alpha (t)h(i^{*},i;t)[{\mathit{X}}^{\mathrm{out}}-{\mathit{\omega }}_i^{\mathrm{out}}]$$

(18)

To sum up, the steps involved in the VQTAM algorithm are as follows:

VQTAM algorithm:

(training part)

Begin

1. Input: X ⁱⁿ, X ^out in training set

2. for $j\leftarrow 0\mathrm{to}\mathrm{length}({\mathit{X}}^{\mathrm{in}})-1$

3. $i^{*}=\underset{i\in A}{\arg min}\{\Vert {{\mathit{X}}_j}^{\mathrm{in}}-{\mathit{\omega }}_i^{\mathrm{in}}\Vert \}$ /* Search activating neuron according to X ⁱⁿ*/

4. ${\mathit{\omega }}_i^{\mathrm{in}}\Leftarrow \alpha (t)h(i^{*},i;t)[{\mathit{X}}^{\mathrm{in}}-{\mathit{\omega }}_i^{\mathrm{in}}]$

5. ${\mathit{\omega }}_i^{\mathrm{out}}\Leftarrow \alpha (t)h(i^{*},i;t)[{\mathit{X}}^{\mathrm{out}}-{\mathit{\omega }}_i^{\mathrm{out}}]$ /* Update the weight vector of neuron */

End

(testing part)

Begin

1. Input: X _test ⁱⁿ in testing set

2. for $j\leftarrow 0\mathrm{to}\mathrm{length}({{\mathit{X}}_{\mathrm{test}}}^{\mathrm{in}})-1$

3. $i\ast =\underset{i\in A}{\arg \min }\left\{\Vert {\left[{\mathit{X}}_{test}^{in}\right]}_j-{\mathit{\omega }}_i^{in}\Vert \right\}$ /* Search activating neuron according to X _test ⁱⁿ*/

4. ${[{\widehat{\mathit{X}}}_{test}^{out}]}_j={\mathit{\omega }}_{i\ast }^{out}$ /* Get the output X _test ^out corresponding to the test data */

5. Output: ${\hat{\mathit{X}}}_{\mathrm{test}}^{\mathrm{out}}$

End

Searching activated neuron by priority search K-means tree algorithm

The weight vector ω _i* ⁱⁿ is the nearest neuron to X ⁱⁿ in the input space, and is searched during the training, testing and estimation of VQTAM, as shown in equation (12). The search algorithm applied has a great impact on the efficiency of VQTAM. The global traversal method needs to calculate all the Euclidean distance from the weight vectors ω _i ⁱⁿ to the input X ⁱⁿ, and searches the minimum value through traversal, which is very inefficient. The K-Dtree algorithm can reduce the time complexity of searching, but for d-dimensional data, the time complexity of K-Dtree data structure search is O(n^d/(d-1)). For high-dimensional data search, especially when the data dimension d>20, K-Dtree data structure cannot provide high efficiency. ω _i ⁱⁿ is 30 dimension. The searching efficiency is severely restricted. The search time complexity of the priority search k-means tree algorithm is O(Ld logn/logK), where n is the amount of data in the data set, and K is the number of nearest neighbors to be searched, and L is the maximum number of data retrieved in the search process. The priority search k-means tree algorithm is suitable for searching high-dimensional data.

The priority search K-means tree algorithm is to divide the space into B different regions using k-means clustering algorithm. Then, the points in each region are partitioned by the same operation until the number of data points in the region is not more than K. In the process of searching, the idea of “divide and conquer” is used to find the nearest point in a smaller area and reduce the amount of data to be searched, so as to improve the efficiency.

The algorithm of establishing a priority search K-means tree data structure is as follows:

K-means tree data structure building algorithm: ⁴²

Begin

1. Input: D← ω _i ⁱⁿ, B, I_max

2. /*B is branch parameter, and I_max is maximum iteration number. */

3. if |D|<B

4.  Create leaf nodes from data sets

5. else

6.  P← uses C_alg algorithm to select B points from data set D

7. for $j\leftarrow 0\mathrm{to}\mathrm{length}(P)$

8.  do

9.   C← clustering data in D centered on P_j

10.   P_j ^new ← finding the mean value of group C data after clustering

11.  while (P = P_j ^new & P_j ^new is the non-leaf node)

12. Output: the entire K-means tree data structure

In the establishment of K-means tree, the algorithm C_alg used to select B points from data set D can be chosen from random selection algorithm, Gonzale algorithm and K-Means++ algorithm. The specific algorithm has little effect on data structure establishment. Generally, random algorithm is recommended. ⁴²

In the built K-means tree, leaf nodes are the data points in the original data set, and non-leaf nodes are the central points of the regions after segmentation. The search process can only be executed in a few areas, thus avoiding the global search of data sets. The K-means tree data structure is searched from the root node. Sub-nodes are sorted according to the distance between the clustering center and the query data points. The nearest sub-nodes are searched in advance. The priority K-means tree search algorithm is as follows:

Priority Search K-means Tree Algorithm ⁴² :

Begin

1. Input: K-means tree, Q, K, L

2. /*Q is the query data X ⁱⁿ ; K is the nearest neighbor number; L is the maximum number of search data */

3. S← the root node /*Initialize a stack S and place the root node in the stack*/

4. do

5. if S(0) is a leaf node

6.  P←S(0) /*Add S(0) to the retrieved data array P */

7. else

8. remove S(0), E←sub-nodes /* The top node goes out of the stack*/

9.  S←push(sort(E))

10. while (isempty(S) & length(P)<L)

11. Output: P(0 to K-1)

End

VQTAM local linear improvement algorithms

The VQTAM algorithm is used to quantify the input space ω ⁱⁿ and the output space ω ^out. The neurons in ω ⁱⁿ and ω ^out correspond through mapping relations. The input X ⁱⁿ is approximated to the nearest neuron ω _i* ⁱⁿ during estimation. The accuracy of the estimation results can be guaranteed when the number of neurons is adequate. However, an increase in the number of neurons increases network size and the decreases the computing efficiency. Thus, local linearization of the activated node is used, which can balance the number of neurons and prediction accuracy. Three improved algorithms for VQTAM are proposed: LLR, LWR and LLE based on local linearization. The inverse kinematics of the robot can be further optimized and ensure the computational efficiency.

Local linearization is performed in VQTAM networks. The K-means tree algorithm is used to search the nearest n data points ω _i*n ⁱⁿ in the input space and the corresponding output data ω _i*n ^out are mapped, as shown in Figure 5.

Figure 5.

Local linearization of VQTAM network.

Improvement algorithm of VQTAM with LLR

It’s assumed that the relationship between local input and output can be expressed linearly, so the LLR can be performed as shown in equation (19).

$$\mathit{YB}=\mathit{Z}$$

(19)

Y is the nearest neighborhood ω _i*n ⁱⁿ of X ⁱⁿ, and Z is the corresponding data point ω _i*n ^out, As shown in equations (20) and (21).

$$\mathit{Y}=\left[\begin{array}{c}{\omega }_{i*1}^{\mathrm{in}}\hfill \\ \multicolumn{1}{c}{⋮}\\ \multicolumn{1}{c}{{\omega }_{i*n}^{\mathrm{in}}}\end{array}\right]\in {\mathbb{R}}^{n\times m}$$

(20)

$$\mathit{Z}=\left[\begin{array}{c}{\omega }_{i*1}^{\mathrm{out}}\hfill \\ \multicolumn{1}{c}{⋮}\\ \multicolumn{1}{c}{{\omega }_{i*n}^{\mathrm{out}}}\end{array}\right]\in {\mathbb{R}}^{n\times l}$$

(21)

The linear relationship between the local input and output can be obtained by solving W . Equation (21) is an over-determined equation, which is generally treated as a least squares problem and solved by a generalized inverse matrix, as shown in equation (22).

$$\mathit{B}=({\mathit{Y}}^T\mathit{Y}{)}^{-1}{\mathit{Y}}^T\mathit{Z}$$

(22)

Y ^T Y is an n-dimensional square matrix. It needs great amount of calculation when the dimension increases. Meanwhile, Y ^T Y may be a singular matrix or ill-conditioned matrix, so it is difficult to find the inverse matrix directly. The SVD (singular value decomposition) can be used to solve the least squares problem. Y can be decomposed as shown in equation (23).

$$\begin{array}{c}\mathrm{SVD}(\mathit{Y})=[\mathit{U}][\mathit{T}][{\mathit{V}}^T]\hfill \\ \mathrm{where}\mathit{Y}\in {\mathbb{R}}^{m\times m},\mathit{U}\in {\mathbb{R}}^{m\times l},\mathit{V}\in {\mathbb{R}}^{l\times l}\hfill \end{array}$$

(23)

The column vectors of U are the left singular vectors of Y , that is, the eigenvectors of YY ^T; the column vectors of V ^T are the right singular vectors of Y , that is, the eigenvectors of Y ^T Y .

$$\mathit{T}=\left[\begin{array}{c}\mathbf{\Sigma }\hfill \\ \multicolumn{1}{c}{\mathit{0}}\end{array}\right]$$

(24)

Where Σ is a diagonal matrix whose value is the singular value υ of matrix Y , that is, the eigenvalue of YY ^T is λ = υ², and U and V ^T are orthogonal matrices.

U can be disassembled as [ U _n , U _m-n ]. For the least squares problem as shown in equation (21), the solution obtained by SVD is shown in Formula (25).

$$\mathit{B}=\mathit{V}{\mathbf{\Sigma }}^{-1}{\mathit{U}}_n^T\mathit{Z}$$

(25)

The input Xⁱⁿ is used to estimate the output, as shown in Formula (26).

$${\hat{\mathit{X}}}^{\mathrm{out}}={\mathit{X}}^{\mathrm{in}}\mathit{B}$$

(26)

Improvement algorithm of VQTAM with LWR

The impact of all neighbor points on the prediction results is consistent in LLR algorithm. Considering the different influence of neighbor points, a distance weight is added in the regression process. A local weighted linear regression (LWR) is performed. The cost function is defined as shown in equation (27).

The weight w_k in equation (27) is determined using the Gauss neighborhood function, as shown in equation (28), where τ is the bandwidth parameter.

$$J(\mathit{\phi })=\sum _{k=1}^n{w}_k{\Vert {\mathit{\omega }}_{i*k}^{\mathrm{out}}-{\mathit{\phi }}_k{\mathit{\omega }}_{i*k}^{\mathrm{in}}\Vert }^2$$

(27)

$$w_k=\exp (-\frac{{\Vert {\mathit{X}}^{\mathrm{in}}-{\mathit{\omega }}_{i*k}^{\mathrm{in}}\Vert }^2}{{\tau }^2})$$

(28)

The analytical solution of the cost function minimization is equivalent to the solution of the over-determined equation shown in equation (29).

$$\mathit{WY}\phi =\mathit{WZ}$$

(29)

W is an n-dimensional diagonal matrix and the diagonal element is w_k.

The SVD for matrix WY is also used to solve the over-determined equation.

$$\mathrm{SVD}(\mathit{WY})=[\mathit{U}][\mathit{T}][{\mathit{V}}^T]$$

(30)

$$\mathit{\phi }=\mathit{V}{\mathbf{\Sigma }}^{-1}{\mathit{U}}_n^T\mathit{WZ}$$

(31)

Improvement algorithm of VQTAM with LLE

It’s assumed that input data can be represented by a linear combination of several samples in its neighborhood. That is to say, the input X ⁱⁿ to be predicted is represented by a linear combination of n data points in its neighborhood in the input space ω _i*n ⁱⁿ .

$${\mathit{X}}^{\mathrm{in}}=\sum _{k=1}^n{c}_k{\mathit{\omega }}_{i*k}^{\mathrm{in}}$$

(32)

The output has the same linear combination.

$${\hat{\mathit{X}}}^{\mathrm{out}}=\sum _{k=1}^n{c}_k{\mathit{\omega }}_{i*k}^{\mathrm{out}}$$

(33)

c_k is the coefficients of the linear combination. The output can be estimated using equation (33) after solving c_k. Define the cost function as shown in equation (34) and rewrite it as a matrix.

$$\begin{array}{c}\multicolumn{1}{c}{J(\mathit{c})=\Vert {\mathit{X}}^{\mathrm{in}}-\sum _{k=1}^n{c}_k{\mathit{\omega }}_{i*k}^{\mathrm{in}}{\Vert }^2}\\ \multicolumn{1}{c}{=\Vert {\mathit{X}}^{\mathrm{in}}-\mathit{Yc}{\Vert }^2}\end{array}$$

(34)

For the least squares problem as shown in equation (33), it can also be solved by SVD of Y .

$$\mathit{c}=\mathit{V}{\mathbf{\Sigma }}^{-1}{\mathit{U}}_n^T{\mathit{X}}^{\mathrm{in}}$$

(35)

Simulation results and discussion

Standard VQTAM network test results

The poses for the sample set are generated according to the robot parameters listed in Table 1, using equations (25) to (30). The generated data can effectively cover the entire workspace of the manipulator.

$$x(t)=1500(1-e^{-\pi t})\cos 1.88\pi t$$

(36)

$$y(t)=1500(1-e^{-\pi t})\sin 1.88\pi t+1500$$

(37)

$$z(t)=2000\cos t+2000$$

(38)

$${\theta }_1(t)=\pi \sin t$$

(39)

$${\theta }_2(t)=2\pi (1-e^{-\pi t})\sin 0.86\pi t$$

(40)

$${\mathit{S}}^1=\left[\begin{array}{ccc}\cos {\theta }_1\sin {\theta }_2\hfill & \sin {\theta }_1\sin {\theta }_2\hfill & \cos {\theta }_2\hfill \end{array}\right]$$

(41)

where t∈[0,1000] and step t = 0.01. Thus Q (t) = [x(t),y(t),z(t); S ¹(t)] is calculated as the output of the sample set. The three components of S ¹ are s_x, s_y and s_z, respectively. p (t) = [l₁(t),l₂(t),l₃(t),l₄(t),l₅(t)] is calculated using the inverse kinematic algorithm. p(t) and Q(t) are used as the sample set. Thus the VQTAM network is trained. The data is selected randomly from the sample set in each epoch. At the same time, in order to test the robustness of the algorithm, random errors are added to the test set data according to N(0,1) distribution.

Table 3 shows the hyper-parameter setting of the VQTAM network training. A VQTAM network of dimension 80×80 (M_x×M_y) is trained based on the hyper-parameters, and the effect of training is evaluated using the test set.

Table 3.

Hyper-parameter setting of VQTAM.

nq	np	epoch	α0	αM	σ0	σM
3	3	5000	0.8	0.001	15	0.001

Root mean squared error (RMSE), R squared (R²), and root mean absolute error (RMAE) are used to evaluate the prediction accuracy of VQTAM for forward kinematics. Table 4 shows the prediction accuracy of the six parameters in output calculated using RMSE, R², and RMAE.

Table 4.

RMSE, R², RMAE of standard VQTAM network.

	x	y	z	sx	sy	sz
RMSE	0.0544	0.0542	0.0689	0.0861	0.0792	0.1425
R2	0.9758	0.9768	0.9600	0.9725	0.9519	0.9649
RMAE	0.7230	0.8522	1.5152	0.7455	1.7239	1.1008

The value of R² is between 0 and 1. The closer it is to 1, the better is the effect of the regression fitting, which means that the learning VQTAM network is more accurate as an approximate model. Generally, R²  > 0.95 can be used in engineering applications. For VQTAM networks of dimension 80×80, the R² values after learning are all above 0.95. This shows that this network is suitable for engineering applications. The RMSE and RMAE reflect the prediction errors in the VQTAM network, and their small values also indicate the prediction accuracy of the network. The convergence result in neural network training process is shown in Figure 6.

Figure 6.

Convergence result. (a) Prediction accuracy of x in the training process. (b) Prediction accuracy of y in the training process. (c) Prediction accuracy of z in the training process. (d) Prediction accuracy of s_x in the training process. (e) Prediction accuracy of s_y in the training process. (f) Prediction accuracy of s_z in the training process.

It can be seen from Figure 6 that the R² of the VQTAM neural network tested using the test set in the training process gradually increases and converges near 1. This confirms the effectiveness of the VQTAM neural network training algorithm.

The consistency between the actual data and the estimated data can be displayed more intuitively through box diagrams, as shown in Figure 7. The box diagrams of the actual values and predicted values are notably similar, which verifies the excellent prediction ability of the VQTAM network of dimension 80×80.

Figure 7.

VQTAM network prediction effect box diagram.

VQTAM local linear improvement algorithms

The output of the inverse kinematics problem is estimated using the VQTAM local linear improvement algorithms. The priority search K-means tree search algorithm is used to search the neighbor data of the input. The influence of parameter k on the prediction accuracy is analyzed, considering a VQTAM network with dimension 70×70 as an example, as shown in Figure 8.

Figure 8.

The influence of the k-nearest neighbor number in the improvement algorithm of VQTAM with LLR.

Figure 8 shows that R² of each parameter increases continuously with the increase in k, at first; and it gradually increases to more than 0.95. This indicates the influence on the prediction accuracy. When k = 7, the overall prediction accuracy is the best.

Figure 9 shows that R² for all variables is more than 0.95 when k = 3. Comparing Figure 8 with Figure 9, the overall prediction accuracy of LWR algorithm is higher than that of LLR algorithm. Comparing among Figures 8 to 10, the overall prediction accuracy of improvement algorithm of VQTAM with LLE is the highest. When k = 7, the overall prediction accuracy is the best. In conclusion, when the number of neurons is small, the prediction accuracy of network can be improved by combining local linear algorithm.

Figure 9.

The influence of the k-nearest neighbor number in improvement algorithm of VQTAM with LWR.

Figure 10.

The influence of k-nearest neighbor number in improvement algorithm of VQTAM with LLE.

Overall examples and test results

According to the VQTAM network parameters shown in Table 3, the prediction accuracy of four algorithms under different network dimensions (30×30, 40×40, 50×50, 60×60, 70×70, 80×80) is studied, as shown in Figure 11. The number of the nearest nodes of the local linear algorithm is set to k = 7.

Figure 11.

Prediction accuracy of 4 algorithms. (a) Prediction accuracy of standard VQTAM network. (b) Prediction accuracy of improvement algorithm of VQTAM with LLR. (c) Prediction accuracy of improvement algorithm of VQTAM with LWR. (d) Prediction accuracy of improvement algorithm of VQTAM with LLE.

Figure 11 shows that an increase of the number of neurons results in a significant increase in the prediction accuracy of the four algorithms. The overall prediction accuracy of improvement algorithm of VQTAM with LLE is the highest. Local linearization of VQTAM yields remarkable results. Using the improved VQTAM algorithm, the prediction accuracy of a network of dimension 50×50 is found to be close to that of the standard VQTAM network of dimension 70×70. Considering improvement algorithm of VQTAM with LLE as an example, the output of the forward kinematics is estimated, as shown in Figure 12. The estimated poses are compared with the actual poses. This shows that the VQTAM neural network can effectively estimate the output

Figure 12.

Forward kinematics output estimation. (a) Forward kinematics output estimation of x. (b) Forward kinematics output estimation of y. (c) Forward kinematics output estimation of z. (d) Forward kinematics output estimation of s_x. (e) Forward kinematics output estimation of s_y. (f) Forward kinematics output estimation of s_z.

The error of VQTAM can be obtained by comparing the result with the actual pose of the robot. The RSME for position estimation is less than 0.1mm, as shown in Table 5. This precision has met the control requirements of general robots in industrial applications.

Table 5.

The RSME of forward kinematics output estimation.

Parameter	x(mm)	y(mm)	z(mm)
RSME	0.0352	0.0327	0.0418
Parameter	sx	sy	sz
RSME	0.0494	0.0443	0.0918

Conclusion

The kinematic of 4SPRR-SPR parallel robot is studied, and the inverse kinematics solution by a geometric method is derived. Given the input pose of the robot, this method can be used to calculate the lengths of the drive limbs conveniently and quickly. It can also be used as the computational basis of VQTAM network to establish sample sets.

A fast forward kinematic mapping method for parallel robots using VQTAM for non-linear dynamic systems identification is proposed. The VQTAM algorithm is constructed to train and test the neuron network. Based on the priority K-means tree search algorithm, the training, testing and estimating process of VQTAM network are improved to enhance the search efficiency of nearest neighbors.

According to local neurons activation, the VQTAM network is improved based on local linearization. To further optimize the prediction accuracy of the network, reduce the dimension of the network in application, the improved local linear algorithms for VQTAM networks are proposed.

The sample data set for VQTAM training is constructed. The VQTAM algorithm and its improved algorithm are tested. The test results show that the increase of network dimension can improve the prediction accuracy of VQTAM, but the computational efficiency is affected. By using VQTAM local linearization improvement algorithms, the estimation accuracy can be optimized for low-dimensional network. The above algorithm can show good prediction accuracy, and the prediction effect of LLE algorithm is the best among the three algorithms.