Optimization Approach for Inverse Kinematic Solution

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.

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,

$$y\left(t\right)=f\left(\theta \left(t\right)\right)$$ E1

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,

$$\theta \left(t\right)=f^{\text{'}}\left(y\left(t\right)\right)$$ E2

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.

$$h=cos\left(\frac{\theta }{2}\right)+i^{\ast }sin\left(\frac{\theta }{2}\right)+j^{\ast }sin\left(\frac{\theta }{2}\right)+k^{\ast }sin\left(\frac{\theta }{2}\right)$$ E3

Figure 1.

Rotation representation of point.

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₁ can be seen rotated through angle θ to point P₃ which crosses the mid-point P₂. Therefore P₁ point is transforming to P₃ 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₁ and h₂ then the composite rotations h₁^∗h₂ can be given by From Figure 1, p¹ is the point vector representing initial position and p³ is the point vector final condition to be transformed. Therefore,

$$\begin{array}{ll}{p}^3& ={h_2}^{\ast }{\left({h_1}^{\ast }{p}^{1\ast }{h}_1^{-1}\right)}^{\ast }{h}_2^{-1}\\ & ={\left({h_2}^{\ast }{h}_1\right)}^{\ast }{p}^{1\ast }\left({h_1^{-1}}^{\ast }{h}_2^{-1}\right)\\ & ={\left({h_2}^{\ast }{h}_1\right)}^{\ast }{p}^{1\ast }\left({h_2^{-1}}^{\ast }{h}_1^{-1}\right)\end{array}$$ E4

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

$$t_r=h+p^1$$ E5

Quaternion transform can be given by,

$$p^2=h^{\ast }{p}^{1\ast }{h}^{-1}$$ E6

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

$$H^{-1}=\left[<h^{-1}>,<-h^{-1}\otimes P\otimes h>\right]$$ E7

where, −h⁻¹ ⊗ P ⊗ h =  − P + [−2k(j × (−P)) + 2j × (j × (−P))] quaternion product, which is defined in the most general form for two quaternions h₁ = (r₁, j₁), and h₂ = (r₂, j₂) as

$$h_1\otimes h_2=\left[r_1{r}_2-j_1•j_2{r}_1{j}_2+r_2{j}_1+j_2\times j_1\right]$$ E8

where j₁ • j₂ and j₂ × j₁ denote dot and cross products, between the 3-dimensional vectors j₁ and j₂. 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

$$H_1\otimes H_2=\left(h_1,P^1\right)\otimes \left(h_2,P^2\right)={h_1}^{\ast }{h}_2,{h_1}^{\ast }{P}^{2 \ast }{h}_1^{-1}+P^1$$ E9

where, ${h_1}^{\ast }{P}^{2 \ast }{h}_1^{-1}=P^2+2r_1\left(j_1\times P^2\right)+2j_1\times \left(j_1\times P^2\right)$

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 θ₁, θ₂, θ₃, θ₄ and θ₅ joint angles for are articulated arm and d₁, d₂ and d₃ are the link offset. a₁, and a₂ represents link lengths.

Figure 2.

(a) Base frame and model of 5R manipulator; (b) configuration of 5R manipulator.

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

$$H_1=\left[<\overline{C_1}+\overline{S_1}\widehat{k}>,<a_1{C}_1\widehat{i}+a_1{S}_1\widehat{j}+d_1\widehat{k}>\right]$$ E10

$$H_2=\left[<\overline{C_2}+\overline{S_2}\widehat{j}>,<-a_2{S}_2\widehat{i}-a_2{C}_2\widehat{k}>\right]$$ E11

$$H_3=\left[<\overline{C_3}+\overline{S_3}\widehat{j}>,<-d_4{S}_3\widehat{i}-d_4{C}_3\widehat{k}>\right]$$ E12

$$H_4=\left[<\overline{C_4}+\overline{S_4}\widehat{i}>,<d_4\widehat{i}>\right]$$ E13

$$H_5=\left[<\overline{C_5}+\overline{S_5}\widehat{j}>,<-d_6{S}_5\widehat{i}-d_6{C}_5\widehat{k}>\right]$$ E14

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

$$H_1^{-1}=\left[<\overline{C_1}-\overline{S_1}\widehat{k}>,<-a_1\widehat{i}>\right]$$ E15

$$H_2^{-1}=\left[<\overline{C_2}-\overline{S_2}\widehat{j}>,<a_2\widehat{k}>\right]$$ E16

$$H_3^{-1}=\left[<\overline{C_3}-\overline{S_3}\widehat{j}>,<d_4\widehat{k}>\right]$$ E17

$$H_4^{-1}=\left[<\overline{C_4}-\overline{S_4}\widehat{i}>,<-d_4\widehat{i}>\right]$$ E18

$$H_5^{-1}=\left[<\overline{C_5}-\overline{S_5}\widehat{j}>,<d_6\widehat{k}>\right]$$ E19

$$Q_i=H_i\otimes H_{i+1}..\dots \dots H_n$$ E20

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

$$\begin{array}{l}{Q}_i=H_i\otimes H_{i+1}..\dots \dots H_n\hfill \\ Q_5=H_5=\left[<\overline{C_5}-\overline{S_5}\widehat{j}>,<-d_6{S}_5\widehat{i}-d_6{C}_5\widehat{k}>\right]\hfill \\ Q_4=H_4\otimes Q_5=\left[<\overline{C_4}+\overline{S_4}\widehat{i}>,<d_4\widehat{i}>\right]\otimes [<\overline{C_5}+\overline{S_5}\widehat{j}>,\hfill \\ <-d_6{S}_5\widehat{i}-d_6{C}_5\widehat{k}>]\end{array}$$ E21

$$\begin{array}{l}{Q}_4=[<\overline{C_4}\overline{C_5}+\overline{S_4}\overline{C_5}\widehat{i}+\overline{C_4}\overline{S_5}\widehat{j}+\overline{S_4}\overline{S_5}\widehat{k}>,\\ <\left(d_4-d_6{S}_5\right)\widehat{i}+d_6{C}_5{S}_4\widehat{j}-d_6{C}_4{C}_5\widehat{k}>]\end{array}$$ E22

$$\begin{array}{l}{Q}_3=H_3\otimes Q_4=\left[<\overline{C_3}+\overline{S_3}\widehat{j}>,<-d_4{S}_3\widehat{i}-d_4{C}_3\widehat{k}>\right]\otimes \\ [<\overline{C_4}\overline{C_5}+\overline{S_4}\overline{C_5}\widehat{i}+\overline{C_4}\overline{S_5}\widehat{j}+\overline{S_4}\overline{S_5}\widehat{k}>,\\ <\left(d_4-d_6{S}_5\right)\widehat{i}+d_6{C}_5{S}_4\widehat{j}-d_6{C}_4{C}_5\widehat{k}>]\end{array}$$ E23

$$\begin{array}{l}{Q}_3=[<\overline{C_4}\overline{C_{3+5}}+\overline{S_4}\overline{C_{3-5}}\widehat{i}+\overline{C_4}\overline{S_{3+5}}\widehat{j}-\overline{S_4}\overline{S_{3+5}}\widehat{k}>,\\ <\left(d_4+d_4-d_6{C}_4{C}_5{S}_3-d_4{C}_3-d_6{S}_5{C}_3-d_4{S}_3\right)\widehat{i}\\ +(d_6{C}_5{S}_4\widehat{j}+(-d_4{S}_3+d_6{S}_5{S}_3-d_6{C}_4{C}_5{C}_3\widehat{-d_4{C}_3)k}>]\end{array}$$ E24

$$\begin{array}{l}{Q}_2=H_2\otimes Q_3=\left[<\overline{C_2}+\overline{S_2}\widehat{j}>,<-a_2{S}_2\widehat{i}-a_2{C}_2\widehat{k}>\right]\otimes \\ [<\overline{C_4}\overline{C_{3+5}}+\overline{S_4}\overline{C_{3-5}}\widehat{i}+\overline{C_4}\overline{S_{3+5}}\widehat{j}-\overline{S_4}\overline{S_{3+5}}\widehat{k}>,\\ <\left(d_4+d_4-d_6{C}_4{C}_5{S}_3-d_4{C}_3-d_6{S}_5{C}_3-d_4{S}_3\right)\widehat{i}\\ +(d_6{C}_5{S}_4\widehat{j}+(-d_4{S}_3+d_6{S}_5{S}_3-d_6{C}_4{C}_5{C}_3\widehat{-d_4{C}_3)k}>]\end{array}$$ E25

Therefore,

$$\begin{array}{l}{Q}_2=\left[<\overline{(C_2}\overline{C_4}\overline{C_{3+5}}-\overline{S_2}\overline{C_4}\overline{S_{3+5}}\right)+\left(\overline{C_2}\overline{S_4}\overline{C_{3-5}}-\overline{S_2}\overline{S_4}\overline{S_{3+5}}\right)\widehat{i}\\ +(\overline{C_2}\overline{C_4}\overline{S_{3+5}}+\overline{S_2}\overline{C_4}\overline{C_{3+5}})\widehat{j}+\left(-\overline{C_2}\overline{S_4}\overline{S_{3+5}}-\overline{S_2}\overline{S_4}\overline{S_{3-5}}\right)\widehat{k}>,\\ <(-a_2{S}_2+d_4{C}_2+d_4{C}_2-d_4{C}_{2-3}+d_6{S}_5{C}_{2+3}\\ -d_4{S}_{2+3}-d_6{C}_4{C}_5{S}_{2+3}-d_4{S}_3)\widehat{i}+(d_6{C}_5{S}_4\widehat{j}\\ +\left(-a_2{C}_2-d_4{S}_2-d_4{S}_2+d_6{C}_4{C}_5{C}_{2+3}+d_4{S}_{2-3}+d_6{S}_5{S}_{2+3}+d_4{\widehat{C}}_{2-3})k>\right]\end{array}$$ E26

$$\begin{array}{l}{Q}_1=H_1\otimes Q_2=\left[<\overline{C_1}+\overline{S_1}\widehat{k}>,<a_1{C}_1\widehat{i}+a_1{S}_1\widehat{j}+d_1\widehat{k}>\right]\otimes \\ \left[<\overline{(C_2}\overline{C_4}\overline{C_{3+5}}-\overline{S_2}\overline{C_4}\overline{S_{3+5}}\right))+\left(\overline{C_2}\overline{S_4}\overline{C_{3-5}}-\overline{S_2}\overline{S_4}\overline{S_{3+5}}\right)\widehat{i}\\ +(\overline{C_2}\overline{C_4}\overline{S_{3+5}}+\overline{S_2}\overline{C_4}\overline{C_{3+5}})\widehat{j}+\left(-\overline{C_2}\overline{S_4}\overline{S_{3+5}}-\overline{S_2}\overline{S_4}\overline{S_{3-5}}\right)\widehat{k}>,\\ <(-a_2{S}_2+d_4{C}_2+d_4{C}_2-d_4{C}_{2-3}+d_6{S}_5{C}_{2+3}-d_4{S}_{2+3}\\ -d_6{C}_4{C}_5{S}_{2+3}-d_4{S}_3)\widehat{i}+(d_6{C}_5{S}_4\widehat{j}+\\ \left.(-a_2{C}_2-d_4{S}_2-d_4{S}_2+d_6{C}_4{C}_5{C}_{2+3}+d_4{S}_{2-3}+d_6{S}_5{S}_{2+3}+d_4{\widehat{C}}_{2-3})k>\right]\end{array}$$ E27

Therefore,

$$\begin{array}{l}{Q}_1=\left[<\overline{(C_1}\overline{C_4}\overline{C_{2+3+5}}+\overline{S_1}\overline{S_4}\overline{S_{2+3-5}})\right)+\left(\overline{C_1}\overline{S_4}\overline{C_{2+3-5}}-\overline{S_1}\overline{C_4}\overline{S_{2+3+5}}\right)\widehat{i}\\ +(\overline{C_1}\overline{C_4}\overline{S_{2+3+5}}+\overline{S_1}\overline{S_4}\overline{C_{2+3-5}})\widehat{j}+\left(\overline{S_1}\overline{C_4}\overline{C_{2+3+5}}-\overline{C_1}\overline{S_4}\overline{S_{2+3-5}}\right)\widehat{k}>,\\ <(a_1{C}_1-a_2{S}_2+d_4{C}_2{C}_1+d_4{C}_2{C}_1-d_4{C}_{2-3}{C}_1+d_6{S}_5{C}_{2+3}{C}_1-d_4{S}_{2+3}{C}_1\\ +d_6{C}_4{C}_5{S}_{2+3}{C}_1-d_6{C}_5{S}_4{S}_1)\widehat{i}+(d_6{C}_5{S}_4+d_6{C}_5{S}_4{C}_1-a_2{S}_2\\ +a_2{S}_2{C}_1-a_2{S}_2{S}_1+d_4{C}_2{S}_1+d_4{C}_2{S}_1-d_4{C}_{2-3}{S}_1+d_6{S}_5{C}_{2+3}{S}_1\\ -d_4{S}_{2+3}{S}_1-d_6{C}_4{C}_5{S}_{2+3}{S}_1)\widehat{j}+(d_1-a_2{C}_2-d_4{S}_2-d_4{S}_2\\ +d_6{C}_4{C}_5{C}_{2+3}+d_4{S}_{2-3}+d_6{S}_5{S}_{2+3}+d_4{\widehat{C}}_{2-3})k>]\end{array}$$ E28

$$O_{j+1}=H_j^{-1}\otimes O_j$$ E29

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

$$\left[R_{\mathit{be}},T_{\mathit{be}}\right]=O_1=\left[<w+a\widehat{i}+b\widehat{j}+c\widehat{k}>,<X\widehat{i}+Y\widehat{j}+Z\widehat{k}>\right]$$ E30

Now using Eq. (29), O₂ will be given by,

$$\begin{array}{l}{O}_2=H_1^{-1}\otimes O_1\hfill \\ \begin{array}{l}{O}_2=\left[<\overline{C_1}-\overline{S_1}\widehat{k}>,<-a_1\widehat{i}>\right]\otimes \\ \left[<w+a\widehat{i}+b\widehat{j}+c\widehat{k}>,<X\widehat{i}+Y\widehat{j}+Z\widehat{k}>\right]\end{array}\hfill \\ \begin{array}{l}{O}_2=[<\left(o_{21}\right)+\left(o_{22}\right)\widehat{i}+\left(o_{23}\right)\widehat{j}+\left(o_{24}\right)\widehat{k}>,\\ <\left(o_{25}\right)\widehat{i}+\left(o_{26}\right)\widehat{j}+o_{27}\widehat{k}>]\end{array}\hfill \end{array}$$ E31

where, $o_{21}=w\overline{C_1}+c\overline{S_1}$ , $o_{22}=a\overline{C_1}+b\overline{S_1}$ , $o_{23}=b\overline{C_1}-a\overline{S_1}$ , $o_{24}=c\overline{C_1}-w\overline{S_1}$ , $o_{25}=X\overline{C_1}-a_1+Y\overline{S_1}$ , $o_{26}=Y\overline{C_1}-X\overline{S_1}$ , o₂₇ = Z.

Now,

$$\begin{array}{l}{\mathrm{O}}_3={\mathrm{H}}_2^{-1}\otimes {\mathrm{O}}_2\hfill \\ {\mathrm{O}}_3=\left[<{\mathrm{o}}_{31}+{\mathrm{o}}_{32}\widehat{i}+{\mathrm{o}}_{33}\widehat{j}+{\mathrm{o}}_{34}\widehat{k}>,<{\mathrm{o}}_{35}\widehat{i}+{\mathrm{o}}_{36}\widehat{j}+{\mathrm{p}}_{\mathrm{z}}\widehat{k}\right]\hfill \end{array}$$ E32

where, $o_{31}=\overline{C_2}{o}_{21}+c\left(\overline{C_2}\overline{S_2}\right)-w\left(\overline{S_2}\overline{S_1}\right)$ , $o_{32}=\overline{C_2}{o}_{22}+\overline{S_2}{o}_{23}$ , $o_{33}=\overline{C_2}{o}_{23}-\overline{S_2}{o}_{22}$ , $o_{34}=\overline{C_2}{o}_{24}-\overline{S_2}{o}_{21}$ , o₃₅ =  − ZS₂ + XC₁C₂ + YS₁C₂ − a₁C₂, o₃₆ = YC₁ − XS₁, o₃₇ = a₂ − ZC₂ + XC₁S₂ + YS₁S₂ − a₁S₂.

$$O_4=\left[<o_{41}+o_{42}\widehat{i}+o_{43}\widehat{j}+o_{44}\widehat{k}>,<o_{45}\widehat{i}+o_{27}\widehat{j}+o_{47}\widehat{k}\right]$$ E33

where,

$$\begin{array}{c}{o}_{41}=\overline{C_{2+3}}{o}_{21}+\overline{S_{2+3}}{o}_{23}\\ o_{42}=\overline{C_{2+3}}{o}_{22}-\overline{S_{2+3}}{o}_{24}\\ o_{43}=\overline{C_{2+3}}{o}_{23}-\overline{S_{2+3}}{o}_{21}\\ o_{44}=\overline{C_{2+3}}{o}_{24}-\overline{S_{2+3}}{o}_{22}\\ o_{45}=-{\mathit{ZS}}_2-a_2{S}_3-{\mathit{XC}}_1{S}_2{S}_3-{\mathit{YS}}_1{S}_2{S}_3+a_1{S}_2{S}_3-\\ {\mathit{ZC}}_2{S}_3+Z-{\mathit{ZC}}_3+{\mathit{XC}}_1{C}_2{C}_3+{\mathit{XS}}_1{C}_2{C}_3-a_1{C}_2{C}_3\\ o_{46}={\mathit{YC}}_1-{\mathit{XS}}_1\\ o_{47}=-{\mathit{ZS}}_2{S}_3+a_2{C}_3+{\mathit{XC}}_1{C}_2{S}_3+{\mathit{YS}}_1{C}_2{S}_3-\\ a_1{C}_2{S}_3+{\mathit{ZC}}_2{C}_3+{\mathit{XC}}_1{S}_2{C}_3+{\mathit{YS}}_1{S}_2{C}_3-a_2{S}_2{C}_3\end{array}$$

Therefore, all the joint variables can be calculated by equating quaternion vector products and quaternion vector pairs i.e. Q₁, Q₂ and Q₃ to O₁, O₂ and O₃ respectively.

$$\Rightarrow {\theta }_1=atan\left[\begin{array}{l}\frac{Y-u_y}{X-u_x}\\ \times \frac{\left(a_1+d_4{C}_2+d_4{C}_2-d_4{C}_{2-3}+d_6{S}_5{C}_{2+3}-d_4{S}_{2+3}+d_6{C}_4{C}_5{S}_{2+3}\right)}{\left(-a_2{S}_2+d_4{C}_2+d_4{C}_2-d_4{C}_{2-3}+d_6{S}_5{C}_{2+3}-d_4{S}_{2+3}-d_6{C}_4{C}_5{S}_{2+3}\right)}\end{array}\right]$$ E34

$${\theta }_2=atan2\left[\left(\frac{Z+v_x}{d_4}\right),\mp \sqrt{1-{\left\{\left(\frac{Z+v_x}{d_4}\right)\right\}}^2}\right]$$ E35

where, d₄S₂ + d₄S₂ − v_x = Z, d₄S₂ + d₄S₂ = Z + v_x, $S_2=\left(\frac{Z+v_x}{d_4}\right)$

$${\theta }_3=atan2\left[\begin{array}{l}\mathrm{m}\sqrt{1-{\left\{\left({\displaystyle \frac{-{\mathit{XS}}_2+{\mathit{XC}}_1{C}_2+{\mathit{YS}}_1{C}_2-a_1{C}_2-v_y}{\left(-d_4-d_6{S}_5\right)}}\right)\right\}}^2},\\ \left(\frac{-{\mathit{ZS}}_2+{\mathit{XC}}_1{C}_2+{\mathit{YS}}_1{C}_2-a_1{C}_2-v_y}{\left(-d_4-d_6{S}_5\right)}\right),\end{array}\right]$$ E36

$${\theta }_4=atan2\left[\begin{array}{l}\mathrm{m}\sqrt{1-{\left\{\left({\displaystyle \frac{-{\mathit{ZS}}_2{S}_3+a_2{C}_3+{\mathit{XC}}_1{C}_2{S}_3+{\mathit{YS}}_1{C}_2{S}_3-a_1{C}_2{S}_3+{\mathit{ZC}}_2{C}_3+{\mathit{XC}}_1{S}_2{C}_3+{\mathit{YS}}_1{S}_2{C}_3-a_2{S}_2{C}_3}{-d_6{C}_5}}\right)\right\}}^2}\\ ,\left(\frac{-{\mathit{ZS}}_2{S}_3+a_2{C}_3+{\mathit{XC}}_1{C}_2{S}_3+{\mathit{YS}}_1{C}_2{S}_3-a_1{C}_2{S}_3+{\mathit{ZC}}_2{C}_3+{\mathit{XC}}_1{S}_2{C}_3+{\mathit{YS}}_1{S}_2{C}_3-a_2{S}_2{C}_3}{-d_6{C}_5}\right),\end{array}\right]$$ E37

$${\theta }_5=atan2\left[\begin{array}{l}\left(\frac{-{\mathit{XS}}_2-a_2{S}_3-{\mathit{XC}}_1{S}_2{S}_3-{\mathit{YS}}_1{S}_2{S}_3+a_1{S}_2{S}_3-{\mathit{ZC}}_2{S}_3+Z-{\mathit{ZC}}_3+{\mathit{XC}}_1{C}_2{C}_3+{\mathit{YS}}_1{C}_2{C}_3-a_1{C}_2{C}_3-d_4}{-d_6}\right),\\ \mathrm{m}\sqrt{1-{\left\{\left({\displaystyle \frac{-{\mathit{ZS}}_2-a_2{S}_3-{\mathit{XC}}_1{S}_2{S}_3-{\mathit{YS}}_1{S}_2{S}_3+a_1{S}_2{S}_3-{\mathit{ZC}}_2{S}_3+Z-{\mathit{ZC}}_3+{\mathit{XC}}_1{C}_2{C}_3+{\mathit{YS}}_1{C}_2{C}_3-a_1{C}_2{C}_3-d_4}{-d_6}}\right)\right\}}^2}\end{array}\right]$$ E38

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 (θ₁, θ₂, θ₃, θ₄ and θ₅) 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):

$$P_c=\left[X_c,Y_c,Z_c\right]$$ E39

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

$$P_d=\left[X_d,Y_d,Z_d\right]$$ E40

Current position of end effector will be compared with the desired position P_d. 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 P_d evaluated by number of iterations.

$$P_{\min }={‖P_d-P_c\left(i\right)‖}^2$$ E41

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

$$P_{\min }=\sqrt{{\left(X_d-X_c\left(i\right)\right)}^2+{\left(Y_d-Y_c\left(i\right)\right)}^2+{\left(Z_d-Z_c\left(i\right)\right)}^2}$$ E42

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.

$${\theta }_c=\left({\theta }_{c1},{\theta }_{c2},{\theta }_{c3},{\theta }_{c4}{\theta }_{c5}\right)$$ E43

$${\theta }_d=\left({\theta }_{d1},{\theta }_{d2},{\theta }_{d3},{\theta }_{d4}{\theta }_{d5}\right)$$ E44

Therefore using square norm the objective function can be given as

$${\theta }_{\min }={‖{\theta }_d-{\theta }_c‖}^2$$ E45

Subjected to joint limits

$${\theta }_1\in \left[{\theta }_{1min},{\theta }_{1max}\right]$$

$${\theta }_2\in \left[{\theta }_{2min},{\theta }_{2max}\right]$$

$${\theta }_3\in \left[{\theta }_{3min},{\theta }_{3max}\right]$$

$${\theta }_4\in \left[{\theta }_{4min},{\theta }_{4max}\right]$$

$${\theta }_5\in \left[{\theta }_{5min},{\theta }_{5max}\right]$$

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

$$\begin{array}{l}F{\left(x\right)}_{\min }=\lambda \left[\sqrt{{\left(X_d-X_c\left(i\right)\right)}^2+{\left(Y_d-Y_c\left(i\right)\right)}^2+{\left(Z_d-Z_c\left(i\right)\right)}^2}\right]\\ +\left[{‖{\theta }_d-{\theta }_c‖}^2\right]\end{array}$$ E46

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₁ = 60 mm, a₂ = 145 mm, d₁ = 150 mm, d₂ = 125 mm, d₃ = 130 mm. Upper and lower limit of five joint angles are: θ₁ = [0, 180]; θ₂ = [0, 150]; θ₃ = [0, 150], θ₄ = [0, 85] and θ₅ = [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

Table 1.

Five different positions and joint variables.

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

Table 2.

Desired joint variables determined through quaternion algebra.

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

Table 3.

TLBO results for joint variable and function value.

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.

Figure 3.

Joint variables and fitness function value for position P1.

Figure 4.

Joint variables and fitness function value for position P2.

Figure 5.

Joint variables and fitness function value for position P3.

Figure 6.

Joint variables and fitness function value for position P4.

Figure 7.

Joint variables and fitness function value for position P5.

Figure 8.

Joint variables and fitness function value for position P1.

Figure 9.

Joint variables and fitness function value for position P2.

Figure 10.

Joint variables and fitness function value for position P3.

Figure 11.

Joint variables and fitness function value for position P4.

Figure 12.

Joint variables and fitness function value for position P5.

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).

Figure 13.

Trajectory of joint angle for quaternion.

Figure 14.

Trajectory of joint angle for TLBO.

Figure 15.

Trajectory of joint angle for GA.

SN	Method	Computational time
1	TLBO	15.671 s
2	GA	7.932 s
3	Quaternion	0.993 s

Table 4.

Computational time for inverse kinematic evaluations.

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.