Five different positions and joint variables.

## Abstract

Inverse kinematics of serial or parallel manipulators can be computed from given Cartesian position and orientation of end effector and reverse of this would yield forward kinematics. Which is nothing but finding out end effector coordinates and angles from given joint angles. Forward kinematics of serial manipulators gives exact solution while inverse kinematics yields number of solutions. The complexity of inverse kinematic solution arises with the increment of degrees of freedom. Therefore it would be desired to adopt optimization techniques. Although the optimization techniques gives number of solution for inverse kinematics problem but it converses the best solution for the minimum function value. The selection of suitable optimization method will provides the global optimization solution, therefore, in this paper proposes quaternion derivation for 5R manipulator inverse kinematic solution which is later compared with teachers learner based optimization (TLBO) and genetic algorithm (GA) for the optimum convergence rate of inverse kinematic solution. An investigation has been made on the accuracies of adopted techniques and total computational time for inverse kinematic evaluations. It is found that TLBO is performing better as compared GA on the basis of fitness function and quaternion algebra gives better computational cost.

### Keywords

- TLBO
- GA
- quaternion
- kinematics

## 1. Introduction

Kinematic chain may consist of rigid/flexible links which are connected with joints or kinematics pair permitting relative motion of the connected bodies. In case of manipulator kinematics it can be categorized into forward and inverse kinematics. Forward kinematics for any serial manipulator is easy and mathematically simple to resolve but in case of inverse kinematics there is no unique solution, generally inverse kinematics gives multiple solutions. Hence, inverse kinematics solution is very much problematic and computationally expensive. For real time control of any configuration manipulator will be expensive and generally it takes long time. Forward kinematics of any manipulator can be understand with translation of position and orientation of end effector from joint space to Cartesian space and opposite of this is known as inverse kinematics. It is essential to calculate preferred joint angles so that the end effector can reach to the desired position and also for designing of the manipulator. Various industrial applications are based on inverse kinematics solutions. In real time environment it is obvious to have joint variables for fast transformation of end effector. For any configuration of industrial robot manipulator for n number of joints the forward kinematics will be given by,

where θ_{i} = θ(t), i = 1, 2, 3, …, n and position variables by y_{j} = y(t), j = 1, 2, 3, …, m.

Inverse kinematics for n number of joints can be computed as,

Inverse kinematics solution of robot manipulators has been considered and developed different solution scheme in last recent year because of their multiple, nonlinear and uncertain solutions. There are different methodologies for solving inverse kinematics for example iterative, algebraic and geometric etc. [8] proposed inverse kinematic solution on the basis of quaternion transformation. [20–36] have proposed application of quaternion algebra for the solution of inverse kinematics problem of different configurations of robot manipulator. [35] presented a quaternion method for the demonstrating kinematics and dynamics of rigid multi-body systems. [34] presented analytical solution of 5-dof manipulator considering singularity analysis. [11] presented quaternion based kinematics and dynamics solution of flexible manipulator. [14] proposed detailed derivation of inverse kinematics using exponential rotational matrices. On the other hand, after numerous surveys on conventional analytical and other Jacobian based inverse kinematics are quite complex as well as computationally exhaustive those are not exactly well suitable for the real time applications. Because of the above-mentioned reasons, various authors adopted optimization based inverse kinematic solution.

Optimization techniques are fruitful for solve inverse kinematics problem for different configurations of manipulator as well as spatial mechanisms. Conventional approaches such as Newton-Raphson can be used for nonlinear kinematic problems and predictor corrector type methods can compute differential problem of manipulator. But major drawback of these methods are Singularity or ill condition which converse to local solutions. Moreover, when initial guessing is not accurate then the method becomes unstable and does not converse to optimum solution. Therefore, recently developed metaheuristic techniques can be used to overcome the conventional optimization drawbacks. Literature survey shows the efficiency of these metaheuristic algorithms or bi-inspired optimization techniques are more convenient to achieve global optimum solutions. The major issue with these nature inspired algorithms is framing of objective function. Even these algorithms are direct search algorithms which do not require any gradient or differentiation of objective function. The comparison of the metaheuristic algorithm with heuristic algorithms is based on the convergence rate as it has been proved that the convergence of heuristic-based techniques is slower. Therefore, to adopt metaheuristic techniques such as GA, BBO, teachers learner based optimization (TLBO), ABC, ACO etc. will be suitable for enhancing the convergence rate and yielding global solution. From literature survey the teaching learning based optimization (TLBO) is similar to swarm based optimization in which the impact of learning methods from teacher to student and student to student has been highlighted. Wherein, the population or swarm is represented by group of students and these students gain knowledge from either teacher or students. If these student gain knowledge from teacher then it is called as teachers phase similarly when students learns form student then it is student phase. The output is considered as result or grades of students. Therefore, number for number of subjects resembles the variables of the function and grades or results gives fitness value, [5, 6]. There are numerous other population centered methods which have been effectively applied and shown efficiency [33]. However, all algorithms are not suitable for complex problem as proved by Wolpert and Macready. On the other hand, evolutionary strategy (ES) based methods such as GA, BBO etc. gives better results for various problems and these methods are also population based metaheuristic [16, 28]. Moreover [22] proposed inverse kinematic solution of redundant manipulator using modified genetic algorithm considering joint displacement (*Δθ*) error minimization and the positional error of end effector. [32] proposed inverse kinematic solution of PUMA 560 robot using cyclic coordinate descent (CCD) and Broyden-Fletcher-Shanno (BFS) technique. [23] proposed IK solution of 4-dof PUMA manipulator using genetic algorithm. This paper uses two different objective functions which are based on end-effector displacement and joint variable rotations. [18] proposed trajectory planning of 3-dof revolute manipulator using evolutionary algorithm. [25] proposed inverse kinematics solution and trajectory planning for D-joint robot manipulator based on deterministic global optimization based method. [1] proposed inverse kinematic solution of redundant manipulator using novel developed global optimization algorithm. [4] proposed inverse kinematic solution of PUMA robot manipulator using genetic programming. In this work, mathematical modeling is evolved using genetic programming through given direct kinematic equations. [17] proposed optimization of design parameter i.e. link length using for 2-dof manipulator. [15] proposed inverse kinematic solution of 2-dof articulated robot manipulator using real coded genetic algorithm. [19] proposed inverse kinematic solution scheme of 3-dof redundant manipulator based on reach hierarchy method. [30] proposed inverse kinematic solution of 3-dof PUMA manipulator for the major displacement propose. In this work they have adopted genetic algorithm with adaptive niching and clustering. [12] proposed inverse kinematic solution of 6-dof MOTOMAN robot manipulator for positioning of the end-effector. In this work they have adopted adaptive genetic algorithm for optimum placement of the end effector. [26] proposed inverse kinematic and trajectory generation of humanoid arm manipulator using forward recursion with backward cycle computation method. [21] proposed inverse kinematic solution for 6R revolute manipulator using real time optimization algorithm. [24] proposed kinematic solution using three different methods such as bee algorithm, neural network which is later optimized by bee algorithm and evolutionary algorithm. [2] proposed kinematic solution of 3-dof serial robot manipulator using real time genetic algorithm. [13] proposed inverse kinematic solution of 6-dof robot manipulator using immune genetic algorithm. [9] proposed conventional approach i.e. penalty function based optimization method for solving IK. Even though few methods can solve hard NP problems, but it requires high-performance computing system and intricate computer programming.

On the other hand, the use of optimization algorithms is not new in the field of multi-objective and NP-hard problem to arrive at a very reasonable optimized solution, the TLBO algorithm have not been tried to solve an inverse kinematics problems and trajectory of joint variables for robot manipulator. Moreover, computational cost for yielding the inverse kinematics solution with adopted algorithms has been compared without any specialized tuning of concern parameters. Therefore, the key purpose of this work is focused on minimizing the Euclidian distance of end effector position based resolution of inverse kinematics problem with comparison of GA and TLBO obtained solution for 5R robot manipulator. The results of all algorithm are computed from inverse kinematics equations and obtained resultant error for data statistics. In other words, end-effector coordinates utilized as an input for joint angle calculations. At the end 4th order spline formula is considered for generation of end effector trajectory and analogous joint angles of robotic arm using TLBO, GA and quaternion. The sectional organization of the paper henceforth is as follows: Section 2 pertains to the mathematical modeling of the 5R robot manipulator and detail derivation of forward and inverse kinematics of 5R manipulator using quaternion algebra. In Section 3 discuss about the inverse kinematic objective function formulation for 5R manipulator. The experimental results as obtained from simulations are discussed elaborately in Section 5.

## 2. Quaternion vector approach for mathematical modeling

This section deals with mathematical modeling of quaternion vector algebra and application for the derivation of inverse kinematics equations. Quaternion vector methods are fruitful for both rotation and translation of a point, line, etc. with references to origin coordinate system irrespective of homogeneous transformation matrix. The Interpolation of series of rotations and translations are quite complex using Euler’s angle method. In other words, the variables lies in isotropic space which is nothing but sphere surface topology and complex in nature. A brief formulation of quaternion mathematics is given in this section for assessment of references and to create background for mathematical derivation of inverse kinematic.

### 2.1. Rotation and translation from quaternion

The above discussions gives the importance of quaternions and the necessity of it. The quaternion rotation and translation are lies in four dimensional space therefore it is quite difficult to represent here or to imagine. Figure 1 represents the rotation through quaternion and Eq. (3) describes rotation of a point in a space mathematically.

The 4-dimensional space, imagination of fourth axis is quite complex. Therefore, in Figure 1 a unit distance point around axis (X, Y, and Z) is given and which traces a circle. When this rotation circle is projected on a plane then the point P_{1} can be seen rotated through angle *θ* to point P_{3} which crosses the mid-point P_{2}. Therefore P_{1} point is transforming to P_{3} following by straight line makes cos(*θ*/2) and sin(*θ*/2).

Now two quaternions can be represented on the basis of above discussed concept. If there is subsequent rotation of two quaternion *h*_{1} and *h*_{2} then the composite rotations *h*_{1}^{∗}*h*_{2} can be given by From Figure 1, *p*^{1} is the point vector representing initial position and *p*^{3} is the point vector final condition to be transformed. Therefore,

Now pure translations *tr* can be done by quaternion operator that is given below,

Quaternion transform can be given by,

Therefore, an expression for the inverse of a quaternion can be given as,

where, −*h*^{−1} ⊗ *P* ⊗ *h* = − *P* + [−2*k*(*j* × (−*P*)) + 2*j* × (*j* × (−*P*))] quaternion product, which is defined in the most general form for two quaternions *h*_{1} = (*r*_{1}, *j*_{1}), and *h*_{2} = (*r*_{2}, *j*_{2}) as

where *j*_{1} • *j*_{2} and *j*_{2} × *j*_{1} denote dot and cross products, between the 3-dimensional vectors j_{1} and j_{2}. Clearly, quaternion multiplication is not commutative. The set of elements {±1, ±i, ±j, ±k} is known as the quaternion group of order 8 in multiplication.

Similarly vector transformation multiplication can be given as

where,

### 2.2. Quaternion derivation for 5R manipulator kinematics

The configuration and base coordinate frame attachment of 5R manipulator is given in Figure 2(a) and MATLAB plot of 5R manipulator is presented in (b). Where *θ*_{1}, *θ*_{2}, *θ*_{3}, *θ*_{4} and *θ*_{5} joint angles for are articulated arm and *d*_{1}, *d*_{2} and *d*_{3} are the link offset. *a*_{1}, and *a*_{2} represents link lengths.

Now quaternion for successive transformation of each joint can be calculated from the Eq. (3) as follows,

Inverse of a dual quaternion can be calculated by Eq. (8),

Where in case of 5R manipulator arm n = 5. Now calculating quaternion vector products using Eq. (20)

Therefore,

Therefore,

Now calculating vector pair of quaternion using Eq. (29), to solve the inverse kinematics problem, the transformation quaternion of end effector of robot manipulator can be defined as

Now using Eq. (29), *O*_{2} will be given by,

where, *o*_{27} = *Z*.

Now,

where, *o*_{35} = − *ZS*_{2} + *XC*_{1}*C*_{2} + *YS*_{1}*C*_{2} − *a*_{1}*C*_{2}, *o*_{36} = *YC*_{1} − *XS*_{1}, *o*_{37} = *a*_{2} − *ZC*_{2} + *XC*_{1}*S*_{2} + *YS*_{1}*S*_{2} − *a*_{1}*S*_{2}.

where,

Therefore, all the joint variables can be calculated by equating quaternion vector products and quaternion vector pairs i.e. *Q*_{1}, *Q*_{2} and *Q*_{3} to *O*_{1}, *O*_{2} and *O*_{3} respectively.

where, *d*_{4}*S*_{2} + *d*_{4}*S*_{2} − *vx* = *Z*, *d*_{4}*S*_{2} + *d*_{4}*S*_{2} = *Z* + *vx*,

## 3. Inverse kinematic solution scheme

In this section optimization algorithms are selected for computation of inverse kinematics solution of 5R manipulator. However, there are various types of optimization algorithms existed and can produce the desired IK solution, the major necessity is to achieve global optimum solution with fast convergence rate. Therefore, selection of appropriate optimization algorithm is important for fitness evaluations and GA is so far best known tool, but on the other hand TLBO has also proven its efficiency and performance. Finally selection of optimization algorithms has been made on global searching point, computational cost and quality of the result.

### 3.1. Optimization approach to solve inverse kinematics

Any Optimization algorithms which are capable of solving various multimodal functions can be implemented to find out the inverse kinematic solutions. The fitness function is given by the Eq. (46) fitness function F(x). Each individual represents a joint variable solution of the inverse kinematic problem for adopted population based metaheuristic algorithm. All individuals moving in D-dimensional search space and sharing the information to find out best fitness value of the function. Each individual contains set of joint angles (*θ*_{1}, *θ*_{2}, *θ*_{3}, *θ*_{4} and *θ*_{5}) of 5R manipulator. The optimum set of joint angle can be find by using appropriate optimization algorithm from given desired position of end effector (X, Y, Z). In case of inverse kinematics of 5R manipulator multiple solutions exist for the single position of the end effector so it is required to find out the best set of joint angle in order to minimize whole movement of manipulator.

For the optimization of joint angle rotation of robot manipulator, one can define objective function or fitness function from joint angle rotation difference and other can be defined from end effector position displacement. These are known as joint angle error and positional function method [3, 7, 10, 29].

#### 3.1.1. Position based function

The current position of the manipulator is described by (39):

Desired position of end effector can be denoted by (40):

Current position of end effector will be compared with the desired position *Pd*. General equation for the fitness function is given in Eq. (41) that is based on the distance norm of homogeneous Euclidian distance between the current positions to the desired position of end effector *Pd* evaluated by number of iterations.

Current position *Pc* can be evaluated from Eqs. (34) through (38). Now putting the value of *Pc* on Eq. (42)

#### 3.1.2. Joint angle error

Corresponding joint error can be given by the difference between current set of joint variables to the final required angles.

Therefore using square norm the objective function can be given as

Subjected to joint limits

Now overall error minimization can be given by using Eqs. (42) and (45),

where *λ* is proportional weight factor for the minimization of the problem and calculation of the entire joint angles base on constraint can be achieved using fitness function (46). The performance of considered algorithm is checked with the parameters: a_{1} = 60 mm, a_{2} = 145 mm, d_{1} = 150 mm, d_{2} = 125 mm, d_{3} = 130 mm. Upper and lower limit of five joint angles are: *θ*_{1} = [0, 180]; *θ*_{2} = [0, 150]; *θ*_{3} = [0, 150], *θ*_{4} = [0, 85] and *θ*_{5} = [15, 45].

## 4. Results and discussions

TLBO and GA has been used to compute the inverse kinematics of 5-R manipulator and comparison of obtained results has been made on the basis of quality and performance. Table 1 gives the five random position of end effector and respective inverse kinematics solutions. Current work is performed in MATLAB R2013a. The data sets are obtained from Eq. (34) through (38). The data sets are generated using quaternion vector based inverse kinematics equations as given in Table 2. These generated data sets are used to compare the IK solution through adopted GA and TLBO. In Table 3, comparative evaluations of fitness function and obtained joint variables through TLBO and GA is presented.

Positions | Joint angles | ||||
---|---|---|---|---|---|

θ_{1} | θ_{2} | θ_{3} | θ_{4} | θ_{5} | |

P1(−76.09, 54.36, −61.94) | 84.559 | 77.518 | 101.74 | 30.616 | 38.697 |

P2(89.69, 192.55, 90.87) | 84.791 | 97.25 | 130.44 | 50.771 | 36.428 |

P3(−4.24, 94.08, 97.55) | 18.384 | 78.688 | 35.234 | 77.708 | 34.889 |

P4(29.10,154.02, −31.52) | 104.43 | 115.47 | 124.11 | 7.3372 | 33.774 |

P5(−184.33, −43.21, 8.27) | 39.177 | 107.13 | 97.052 | 65.672 | 15.374 |

SN | Position of joints determined through quaternion algebra | |||||||
---|---|---|---|---|---|---|---|---|

θ_{1} | θ_{2} | θ_{3} | θ_{4} | θ_{5} | X | Y | Z | |

1 | 112.5641 | 47.3165 | 8.2447 | 65.8373 | 39.8977 | −186.6903 | 183.0670 | −14.7039 |

2 | 153.1316 | 21.9812 | 126.9031 | 57.2629 | 30.4168 | −92.6981 | 32.1423 | 157.3316 |

3 | 66.1779 | 143.2985 | 14.6124 | 73.6231 | 41.5228 | −131.5420 | −22.3866 | −32.2155 |

4 | 57.6085 | 119.6396 | 104.5818 | 71.1946 | 33.8225 | −10.7684 | 111.7435 | 77.4862 |

5 | 31.4308 | 2.9242 | 71.5757 | 63.0358 | 39.8749 | 64.7966 | 172.0372 | 151.5714 |

6 | 124.3702 | 116.7337 | 102.4999 | 53.1482 | 22.4807 | −111.8590 | −59.6708 | 60.8590 |

7 | 89.1765 | 13.1827 | 101.6747 | 80.3340 | 29.4704 | −76.9533 | 96.2813 | 121.3505 |

8 | 5.6698 | 30.9685 | 57.2308 | 29.3079 | 29.6421 | 174.3873 | 107.6283 | 143.1839 |

9 | 131.5857 | 108.7086 | 92.8278 | 5.6664 | 36.4826 | −104.6410 | 109.7511 | 40.5523 |

10 | 32.8579 | 102.3539 | 138.8770 | 26.4141 | 33.7466 | 146.4984 | 48.7416 | 54.4041 |

11 | 134.1878 | 70.4224 | 26.3511 | 82.2471 | 44.0566 | −188.7864 | 15.4108 | −53.9823 |

Positions | TLBO joint angles | Function value | ||||
---|---|---|---|---|---|---|

θ_{1} | θ_{2} | θ_{3} | θ_{4} | θ_{5} | ||

P1(−76.09, 54.36, −61.94) | 86.598 | 72.165 | 72.165 | 40.894 | 30.459 | 0 |

P2(89.69, 192.55, 90.87) | 83.874 | 69.895 | 69.895 | 39.607 | 30.76 | 0 |

P3(−4.24, 94.08, 97.55) | 84.512 | 70.427 | 70.427 | 39.909 | 30.686 | 0 |

P4(29.10,154.02, −31.52) | 85.566 | 71.305 | 71.305 | 40.406 | 30.53 | 0 |

P5(−184.33, −43.21, 8.27) | 87.818 | 73.181 | 73.181 | 41.469 | 30.364 | 0 |

Positions | GA Joint angles | Function value | ||||

θ1 | θ2 | θ3 | θ4 | θ5 | ||

P1(−76.09, 54.36, −61.94) | 60.619 | 49.504 | 58.384 | 62.281 | 27.903 | 0 |

P2(89.69, 192.55, 90.87) | 88.293 | 34.091 | 14.439 | 15.241 | 51.738 | 0 |

P3(−4.24, 94.08, 97.55) | 55.004 | 49.274 | 63.942 | 47.842 | 33.633 | 0 |

P4(29.10,154.02, −31.52) | 72.594 | 22.689 | 68.297 | 85.886 | 27.044 | 0.0137 |

P5(−184.33, −43.21, 8.27) | 25.669 | 70.588 | 31.341 | 66.807 | 52.884 | 0 |

This work does not use special tuning of various parameters of GA and TLBO algorithm. In future research the sensitivity analysis can be performed to achieve better results. From Table 3, TLBO generated solutions for the position 4 is better as compared to GA in account of fitness function evaluation. There are different distance based norms, one of them is Euclidean distance norm and which is used here for minimum distance between the end effector positions. If the distance between two points reached to 0 or less than 0.001 than the evolutions of fitness function can be reached best or global minimum value. It is clear that the obtained fitness value is less than the defined distance norm so adopting these algorithms are fruitful and qualitative.

Figures 3–7 signify the best fitness function value and analogous joint variables for position 1. These figures show efficiency of adopted algorithms for IK solution of 5-R manipulator. The convergence of objective function evaluation lies to zero error for GA and TLBO algorithms while for position 4, GA yields 0.013 error. It means that GA is less performing as compared to TLBO. From Figures 8–12, the results obtained through GA shows in terms of convergence and histogram graph and the obtained joint angles are in radian which is later converted into degree and given in Table 3. The GA results are obtained through MATLAB toolbox and that shows the zero convergence in single run. Figures 8–12, it can be seen that the generated solutions for joint angles are multiple for single position and similarly there are multiple fitness function evaluations. The best fitness function achieved here using the termination criteria and the corresponding joint variables has taken for comparison.

The proposed work is performed in dual core system with 4 GB RAM computer. It has been observe that the convergence of the solution for GA is taking less computation time as compared to TLBO and quaternion algebra. Corresponding joint angles trajectory using 4th order cubic spline is presented in Figure 13. Using inverse kinematic solution joint variables are used to calculate the joint space trajectory for TLBO and GA as presented in Figures 14 and 15. Final time has been taken *t*_{f = 6} second for trajectory generation but to complete this trajectory overall computational time is 5.674 seconds. The computation time for TLBO is 15.671 seconds which is more than the GA i.e. 7.932 seconds. Therefore, on the basis of computational cost GA is performing better than TLBO while quaternion algebra taken least time (Table 4).

SN | Method | Computational time |
---|---|---|

1 | TLBO | 15.671 s |

2 | GA | 7.932 s |

3 | Quaternion | 0.993 s |

## 5. Conclusions

In this paper, the work discourses the problem associated to the optimization of positional and angular error of end effector using TLBO and GA for 5R robot manipulator. Metaheuristic algorithms like PSO, GA, ABC, etc. have been used in various field of industrial robotics but the most critical issue is to solve inverse kinematic problem for any configuration of robot manipulators. Most of the optimization approach are being used for numerical solution but it has been observed that the numerical solutions does not yield solution when the manipulator is in ill-conditioned besides this it has also been observed that classical optimization methods converge in local minima. Therefore in this work global optimization method like TLBO and GA is adopted and after analyzing the results it can be concluded that adopted optimization algorithms convergence rate is higher and complexity does not increase with the manipulator configuration. Although many researchers are tried to obtain global solution but the computations cost are more in the problem henceforth overcoming the problem of computational cost with quaternion objective function.

The adopted algorithms are very much appropriate for constrained and unconstrained problems. To estimate the effectiveness of considered algorithms, comparison has been made with quaternion algebra. Table 3 gives the comparative results of adopted algorithm and proposed quaternion solutions of 5-R manipulator. This work considered forward and inverse kinematic equations for preparing the objective function for TLBO and GA. These adopted algorithms has shown the potential of getting faster convergence and yielding global optimum solution for the stated problem. In future the tuning of various parameters of GA and TLBO can be considered so as to avoid trapping in local minimum point. Even the hybridization of these algorithms may be proposed and adopt for the IK problems.