Use Improved Differential Evolution Algorithms to Handle the Inverse Kinetics Problem for Robots with Residual Degrees of Freedom

In this study, the Self-adaptive strategy algorithm for controlling parameters in Differential Evolution algorithm (ISADE) improved from the Differential Evolution (DE) algorithm, as well as the upgraded version of the algorithms has been applied to solve the Inverse Kinetics (IK) problem for the redundant robot with 7 Degree of Freedom (DoF). The results were compared with 4 other algorithms of DE and Particle Swarm Optimization (PSO) as well as Pro-DE and Pro-PSO algorithms. These algorithms are tested in three different Scenarios for the motion trajectory of the end effector of in the workspace. In the first scenario, the IK results for a single point were obtained. 100 points randomly generated in the robot ’ s workspace was input parameters for Scenario 2, while Scenario 3 used 100 points located on a spline in the robot workspace. The algorithms were compared with each other based on the following criteria: execution time, endpoint distance error, number of generations required and especially quality of the joints ’ variable found. The comparison results showed 2 main points: firstly, the ISADE algorithm gave much better results than the other DE and PSO algorithms based on the criteria of execution time, endpoint accuracy and generation number required. The second point is that when applying Pro-ISADE, Pro-DE and Pro-PSO algorithms, in addition to the ability to significantly improve the above parameters compared to the ISADE, DE and PSO algorithms, it also ensures the quality of solved joints ’ values.


Introduction
The robot Inverse Kinematics problem involves finding the joints' variable values that match input parameters of position and direction of the end effector [1]. These matched variable values will ensure that subsequent robot control will follow the desired trajectory. This is one of the important issues in the robotic field because it is related to other aspects such as motion planning, dynamic analysis and control [2]. Traditionally, there are several methods to resolve inverse kinematics problem for robots such as: geometry method is the method using geometric and trigonometric relationships to solve; the iterative method is often required inversion of a Jacobian matrix, etc. However, when applying these methods to solve the IK problem for robots, especially with redundant robots, it is often much more complicated and time-consuming. The reason is the nonlinearity of the formulas and the geometry between the workspace and the joint space. In addition, the difficult point is in the singularity, the multiple solutions of these formulas as well as the necessary variation of the formulas corresponding to the changes of different robot structures [3][4][5].
In addition to those existing methods of solving the IK problems, in recent years, the application of meta-heuristic optimization algorithms has become increasingly common. 8 optimization algorithms applied in [5] in the cases of a single point or a whole trajectory endpoint. The simulation results showed that the PSO algorithm can effectively solve the IK problem. In [6] the authors used algorithms such as ABC, PSO, and FA to solve the inverse kinematic requirement for Kawasaki RS06L 6-DoF robot in the task of picking and place objects. Ayyıldız et al. compared the results of all IK tests for a 4-DOF serial robot using 4 different algorithms: PSO, QPSO, GA and GSA [7]. Two versions of the PSO algorithm have been used to solve the IK problem for robots with a number of degrees of freedom from 9 to 180 [8]. In recent research [9], Malek et al. used PSO algorithm to handle inverse kinematics for a 7-DoF robot arm manipulator. The study mentioned both the requirements for the location and the direction of the endpoint, however, it only solved for 2 different end effector positions. Laura et al., in [10] used DE algorithm for the IK problem of 7-DoF robot. The problem was solved for specific points, but the quality evaluation parameters such as endpoint position deviation, execution time as well as the values of the joints' variable did not reach impressive quality. Ahmed El-Sherbiny et al. [11] proposed to use ABC variant algorithm for solving inverse kinematics problem in 5 DoFs robot arm. Serkan Dereli et al. [12] used a quantum behave partial algorithm (QPSO) for a 7-DoF serial manipulator and compare the results with other techniques such as firefly algorithm (FA), PSO and ABC. In this study, the self-adaptive control parameters in Differential Evolution (ISADE) algorithm, that developed [13,14] by authors, was applied to solve the problem of inverse kinematic for a 7-DOF serial robot. To compare the results, this IK problem was also handled by applying DE and PSO algorithms. In addition, the study also compared the results in the application of the above algorithms with the search space improvement of joints' variables (Pro-ISADE, Pro-PSO and Pro-DE) [15].
The remainder of the paper is divided into the following sections: Section II describes the experimental model. The theory of the PSO, DE and ISADE algorithms as well as the algorithms with improved search area, Pro-PSO, Pro-DE and Pro-ISADE, will then be presented in Section III. Section IV covers scenarios and object functions that will be applied to calculate the IK. The results after applying the algorithm are shown and compared in Section V. Finally, the conclusions are outlined in Section VI.

Testing model
The residual driven robots have many advantages such as easy escape from obstacles, flexible movement as well as a large operation space. However, their disadvantage is the complexity of the robot structure [16]. In this study, a serial redundant manipulator robot was used to evaluate the algorithm in resolving the inverse kinematics requirements. The simplified robot model was shown in the Figure 1. As in the figure, this serial robot manipulator is of type 7R (R: Revolute). The parameters of the D-H table of the robot are given in Table 1.
The homogeneous transformation matrix can be used to obtain the forward kinematics of the robot manipulator, using the DH parameters in Eq. (1) [17].
where S and C denote the sine and cosine functions. The position and orientation of the end-effector can be determined by Eq. (2): n y s y a y y 5 n z s z a z z 5 0 0 0 1 4 -π/2<q 4 <π/2 0 l 4 = 300 π/2 5 -π/2<q 5 <π/2 0 l 5 = 200 -π/2 6 -π/2<q 6 <π/2 0 l 6 = 200 0 7 -π/2 <q 7 <π/2 d 7 =5 l 7 = 100 0 With: When solving the problem of inverse kinematics, with the endpoint coordinates as on the left side of Eq. (10), we need to find the values of the matching variable q. However, according to the Equation, the number of equations is much less than the number of variables. This makes it very difficult to find a unique and exact matching solution. In this study, ISADE algorithm and ISADE with searching space improvement algorithm (Pro-ISADE) were used to solve the IK problem for the robot. To compare the results, the study also used some other optimization algorithms such as PSO, DE as well as Pro-PSO and Pro-DE to solve the same IK problem for the robot above.

PSO
Particle swarm optimization was developed flying Kenney and Eberhart [18,19] based on observing the moving characteristics of bird flock and fish school. In this algorithm the individual of the population is called particle. The particle of the population (Called swarm) can move in its space and offer a potential solution. Particles can memorize best condition and find and exchange information to other members. Each particle in the population has two characteristics: position and velocity. Starting with the particle population, each particle monitors its coordinates and updates position and speed according to the best solution for each iteration. The velocity and position values are shown in the following equation: In particular, x xi , v i are the position and velocity of the particle i-th, respectively; d is number of dimension; w is the inertia weight factor, c 1 and c 2 are cognitive learning rate and social learning rate, respectively; pi is the pbest value of i_th particle; What is gbest value of the population.

DE
Differential Evolution (DE) algorithm is a population-based stochastic optimization algorithm recently introduced. DE works with two populations; old generation and new generation of the same population. The population is randomly initialized within the initial parameter bounds individuals in the population has two characteristics: position and velocity. Starting with the individual population, each individual monitors its coordinates and updates position and speed according to the best solution for each iteration. Velocity values (V) is randomly created in one of eight ways: In particular, F is Scaling factor, r 1 , r 2 , r 3 , r 4 , r 5 is random solution, r 1 , r 2 , r 3 , r 4 , r 5 ∈ 1, 2, 3, … , N p È É and r 1 6 ¼ r 2 6 ¼ r 3 6 ¼ r 4 6 ¼ r 5 6 ¼ i, X best is population filled with the best member.
Position values new (U) shown in the following Equation: U ¼ X: * FM mpo þ V: * FM mui (13) FM mui are all random numbers <0.9, FM mpo is inverse mask to FM mui:

ISADE
In the [13,14], we suggested to develop a new version of DE algorithm that can automatically adapt the learning strategies and the parameters settings during evolution. The main ideas of the ISADE algorithm are summarized below.

Mutation operator
ISADE probabilistically selects one out of several available learning strategies in the mutation operator for each individual in the current population. In this research, we select three learning strategies in the mutation operator as candidates: "DE/best/1/bin", "DE/best/2/bin" and "DE/rand to best/1/bin" that are respectively expressed as: DE=best=2 DE=rand to best=1 : gare current population and design variable, respectively.
}DE=Randtobest=1=bin} strategy usually demonstrates good diversity while the }DE=best=1=bin} and }DE=best=2=bin}strategy show good convergence property, which we also observe in our trial experiments.

Adaptive scaling factor F and crossover control parameter CR
In the ISADE algorithm, the author suggested to use the sigmoid function to control neighborhood parameter. we sort the particles by estimating their fitness. A ranked particle is labeled with ranked number and assigned F that corresponds with its number. The formula for F by sigmoid function as following: Where: α, idenote the gain of the sigmoid function, particle of the i th in NP, respectively.
For better performance of ISADE, the scale factor F should be high in the beginning to have much exploration and after curtain generation F needs to be small for proper exploitation. Thus, we proposed to calculate the F as follow: Where: F max , F min , iter, iter max and n iter are the lower boundary condition of F, upper boundary condition of F, current generation, maximum generation and nonlinear modulation index, respectively.
The author introduced a novel approach of scale factorF i of each particle with their fitness in Eq. (15). Thus, in one generation the value of F iter Þare not the same for all particles in the population rather they are changed in each generation. The final value of scale factor for each generation is calculated as follow: The control parameter CR is adapted as following: The ISADE algorithm was summarized as in the Figure 2.

Cost functions and Algorithms with searching space improvement
As mention in the introduction part, the disadvantage of many studies using optimization algorithms to solve the IK problem of redundant robots is to focus on the results related to the optimal running process such as execution time, number of generation … but have not yet considered the feasibility of the joints' variable values. In order to overcome these drawbacks, the author of this research [15] proposed this algorithm that is explained as following: The solution to improve the continuity of joints' values constrains the initialization domain of X. This help the program to achieve the dual goal of increasing calculation speed, accuracy and ensuring continuity for the value of joints' variables. In this algorithm, firstly the robot from any position moves to the first point of the trajectory. With this first point, the initialization values for the particles are randomly selected in the full Range of Motion (RoM) of joints. In addition, the target function in this case has the form: where the values q k i and q k i (i = 1) are the joints' variable values at the original position and 1st point on the trajectory, respectively; (x i , y i , z i ) and (x ei , y ei , z ei ) are the End-effector coordinates for the i-point (i = 1) found by the algorithm and the desired End-effector coordinates; (Rx i , Ry i , Rz i ) and (Rx ei , Ry ei , Rz ei ) are corresponding rotation cosine angles performing orientation of the end-effector which are found by Algorithm and orientation of the desired end-effector; a, b are penalty coefficients. Cost function as Eq. (21) ensures the energy spent in the joints to reach the 1st desired position is minimized. Besides, it also minimizes the distance error between the actual and desired end-effector position. The condition to stop for points of trajectory is that the Cost Func.1 value is less than value of e or the number of iterations reaches 600 and the number of times algorithm running <10. After calculating for 1st point of the trajectory, the remaining points are calculated with a search limitation around the previous optimal joints' values. By using this suggested range, the program's search space will be limited while ensuring the continuity of the joint variables. In this case, the target function is still the same as the function of 1st point, but it has coefficient a = 0.

Scenario 1
In Scenario 1, an endpoint in the workspace were randomly selected; the PSO; DE and ISADE algorithms were then applied to solve the required problem. The purpose of this Scenario is to compare the convergence speed of the three algorithms. In this case, since the initial and the desired endpoints can be far apart, the Pro-PSO; Pro-DE and Pro-ISADE algorithms cannot be applied.

Scenario 2
In this case, the robot was required to move the endpoint through 100 points in the robot's working space one after another. These points were selected at random for the purpose of testing the effectiveness of each algorithm with many distinct points. Similar to the previous case, in this Scenario we also only applied the algorithms PSO, DE and ISADE with the solution space of the matching variable which limits the motion of these joints.

Scenario 3
The manipulator robot was required to move the end effector following a certain trajectory. The selected trajectory is spiral, and it is described by the following function: Where:

Experimental setup
The main task of this study is to find the optimal value of the joints' variable to ensure the end effector of robots can reach the desired points. The desired point positions of the Scenario 2 and 3 are shown as the Figure 3. Research using the ISADE and Pro-ISADE algorithm, which were developed by the authors [13][14][15], to get simulation results of inverse kinematics problem and then compared it with the results when using PSO, DE and Pro-PSO, Pro-DE algorithms. When solving the IK problem for the 7-DoF serial robot manipulator, the study focused on three main aspects. The first of these is the sensitivity of the solution -in the other word, the amount distance error of end effector is minimum. The second criterion was the execution time. In order to avoid the endless loop, the maximum numbers of generation iter max ð Þwere set as 600, 600 and 130 for PSO (Pro-PSO), DE (Pro-DE) and ISADE (Pro-ISADE), respectively. And the final aspect is the searching space of joints' variables. Normally, Normally, almost all studies have been using the Range of Motion (RoM) of joints for its boundary space. Our algorithm [15] proposed to use the searching space of current generation is around previous optimal joints' values. In the Table 2, the ubs iþ1 and lbs iþ1 are the joints' upper and lower boundary of the current generation.C 1 and C 2 are weights of personal best and global best, respectively. w is the inertia weight. ρ is the number of run for each algorithm to choose the best result. Besides, after some trial runs for the algorithms, we noticed that our ISDE algorithm gave much better results than DE and the least was the PSO algorithm. Thus, when setting up the maximum distance error by the fitness value setting for the end effector position, the study set the value of 1e À 14 m ð Þ; 1e À 15 m ð Þ and 1e À 17 m ð Þ for PSO (Pro-PSO); DE (Pro-DE) and ISADE (Pro-ISADE), respectively or that can be seen in the Table 2. In this research, the proposed and other methods were tested in the two different Scenarios. Both the first and second Scenario was coded by Matlab version 2019a and run on the computer equipped with an Intel Core i5-4258U @2.4GHz processor and 8 GB Ram memory.

Scenario 1 results
After applying the inverse kinematic problem processing algorithms for a single endpoint, the results are shown in Figures 4 and 5. All algorithms are able to handle  Table 2.

Optimization parameters used in PSO, Pro-PSO, DE Pro-DE, and ISADE, Pro-ISADE.
the inverse kinetics problem, but the best results have been obtained with the ISADE algorithm as shown in Table 3. Figures 4 and 5 show convergence speed of algorithms corresponding to the number of iterations and processing time, respectively. The results show that the   processing speed of the ISADE algorithm is the best, followed by the DE algorithm and finally with the PSO algorithm. In Table 3 the study of selecting stop conditions for algorithms is the maximum number of iterations of 85 rounds. After 10 runs, the best results are shown in the table. The ISADE algorithm gives the best processing results in terms of both quality and speed. The endpoint deviation can reach 2.8422e-13 (m) in 0.049 (s) time. For the PSO algorithm, it can handle the reverse kinematic problem for the end point with an accuracy of 2.6815e-4 in a period of 0.0941 (s). and, 5.7514e-10 (m) and 0.0715 (s) are the accuracy of end effector and execution time for DE algorithm.

Scenario 2 results
As mentioned above, in this Scenario 2, algorithms was used to resolve inverse kinematics problem for 100 randomly chosen points within the workspace of the robot. When processed at each point, the end effector started at the same initial position of [0 0 0 0 0 0 0] for 7 serial joints values. Because the end effector points all come from the same starting point to go to each of the 100 points, the study only used the ISADE algorithm and compares with the results from PSO and DE algorithms without using the Pro-ISADE algorithm as well as Pro-DE and Pro-PSO.
The 100 randomly selected points were shown in the Figure 3a. Results when applying ISADE and the other algorithm were presented in the Figure 6. As shown in the Figure, all algorithms have solved problem well. In particular, with the ISADE algorithm, although the fitness value in experimental setup required 1000 and 100 times higher than the required by applying the PSO and DE algorithms, respectively, it was not only guaranteed required precision but also showed faster processing speed and fewer iterations compared to the 2 other algorithms. Specifically, as shown in Figure 6b and c and especially Table 4, the average execution time when using ISADE to solve IK of each points was around 0.0685 second, while this value of the PSO and DE algorithm were on average 0.2307 (s) and 0.0978 (s) respectively. The main reason for this, as seen in Figure 6b and Table 4, was the number of generations to reach the optimal values much higher in PSO algorithm and slightly higher in DE algorithm, compared to in ISADE algorithm. Specifically, the PSO algorithm needed an average of 413.24 and the DE algorithm needed average of 124.45 loops to find a solution, while the ISADE algorithm used an average of 85.63 loops. Another remarkable thing is although there was not much difference in the number of iterations to solve the problem between the two algorithms DE and ISADE, but the ISADE algorithm still gave a processing speed of 1.42 times higher than DE algorithm though required 100 times more accuracy for the ISADE algorithm. This demonstrated the very high efficiency of the ISADE algorithm when it was applied to handle inverse kinematics problem for this robot. In short, in the optimization study for randomly chosen points in working space, the ISADE algorithm presented the best algorithm to resolve the IK requirement in term of accuracy, iteration and execution time.

Scenario 3 results
In Scenario 2, the end effector moved through the 100 points located on a specific trajectory that was defined in Eq. (22) and shown in Figure 3b. The main difference between Scenario 2 and Scenario 3 is that, instead of after solving each IK problem for each point, the end effector goes back to the original point to continue processing for the next points like in Scenario 2, in Scenario 3 the end effector starts from previous point in order to calculate for the next point. Stemming from this feature, the searching space of joints' variable also starts previous optimal joints' values. However, depending on the searching space we have 2 smaller cases such as: • Scenario 3.1: Searching spaces for joints' variables are RoMs. Then, like the Scenario 2, the study compared the results when using the ISADE algorithm with the results when using the PSO and DE algorithms. • Scenario 3.2: Searching spaces for joints' variables are around the previous optimal joints' values. The study compared the results when using Pro-ISADE algorithm with when using Pro-PSO and Pro-DE algorithms.
The results were presented in the Figure 7 and Table 5. Similar to the Scenario 2, although the experimental installation required the ISADE (and Pro-ISADE) algorithm to be 100 and 1000 times more accurate than the algorithm DE (Pro-DE) and PSO (Pro PSSO), respectively, all of 6 algorithms gave appropriated solutions for all the points in the trajectory. It can be seen that, in both cases 3.1 and 3.2 the ISADE and Pro-ISADE algorithms showed the best ability to resolve the inverse kinematics problems in all 3 aspects: accuracy, execution time and number of generations. More specifically, in Scenario 3.1, when searching space for joints' variables were RoMs, the average achieved accuracies for ISADE was around 2.0748e-14 (m) that is much better than the values of 7.5404e-13 (m) and 2.2260e-13 (m) corresponding for PSO and DE algorithms. Although the ISADE algorithm was set to a fitness value to achieve such higher accuracy, the execution time of the algorithm was still below the time of PSO and DE algorithm. These average execution time values were 0.0679 (s); 0.0845 (s) and 0.3478 (s) second for ISADE, DE and Pro algorithm, respectively. The above results can be partly explained based on the number of necessary iterations that each algorithm was needed to find the optimal values of joints variables. From Figure 6c, it showed that, when solving the IK problem for almost points in the spiral trajectory, the ISADE method used the least number of iterations. The Table 5 presented more clearly, on the average each point in the trajectory the ISADE needed 85.19 generations to find the optimal values, these means number for DE and PSO algorithm are 125.44 and 391.1 In Scenario 3.2, the searching space for joints' variables were around previous optimal values that were set up as in the Table 2. Similar to the Scenario 3.1, all of the comparison parameters gotten from using Pro-ISADE algorithm were better than that values from Pro-DE and Pro-PSO algorithms. These parameters are described in the as well as Table 5. In order to comparison between Scenario 3.1 with Scenario 3.2, all average parameters was shown in the Table 5. From all comparison, the proposed ISADE or Pro-ISADE were always proved the best solution to solve the inverse kinematics requirements for the manipulator robot. Moreover, Table 5 also showed that, the Pro-ISADE had better performance compared to ISADE. By using Pro-ISADE algorithm, it reduced all of parameters including distance error, execution time and number of generations.
Another very important result gotten from Scenario 3.2 is the quality of joints' values. Figure 8 show the joints' value in two cases of using ISADE in Scenario 3.1 and using Pro-ISADE in Scenario 3.2. It is clear that the joints' value in the Scenario 3.1 were change dramatically. On the contrary, the values of joints in Scenario 3.2 changed continuously and slowly. The quality of joints variable values as Figure 9b, that received by using Pro-ISADE, will ensure feasibility in the next stages of calculation and design for the robot. These values, along with the values of speed, acceleration, as well as the weight parameters of the stages, will be used in the dynamic problem as well as in future control.  In short, after comparing the results of Scenario 3.1 and 3.2, it is possible to conclude that the ISADE algorithm and Pro-ISADE are the best solutions to solve the IK problem for the robot in all aspects: endpoint accuracy, execution time and number of generation. The Pro-ISADE algorithm not only guarantees the above parameters, it also ensures the quality of the joints' variables to serve the next computational and design stages. Table 5 summarizes results of the average error, the standard deviation of error (STD), the average iteration and the average execution time of all Scenarios. As in the table, the algorithms of ISADE and Pro-ISADE got the better results than the other algorithms.
As mentioned at the beginning of this article, intelligent optimization techniques have been using more and more popular in difficult and complex tasks including the IK problem for redundant manipulator robots. Table 6 shows some studies used meta-heuristic optimization algorithms to resolve the inverse kinematics task for different robot models. The Table presents: the used algorithm for the IK calculation, selected manipulators for the test, the algorithms that are used to comparison. For example, El-Sherbiny et al. [11] used the Adaptive Neuro Fuzzy Inference System (ANFIS) algorithm to calculate the IK problem of a 5 DOF robot, and then compared results with GA algorithm. Both algorithms could get the appropriate solutions, but ANFIS algorithm proved to be the best one. The comparison also shows that a number of studies [12,22,23], using optimal algorithms such as PSO, ABC, Q-PSO … handle the inverse kinetic requirements for the model of 7 degrees of freedom. All the used algorithms have proven the ability to handle the problem, but it is not difficult to see that most of these studies have the lower accuracy and processing speed than the ISADE as well as the Pro-ISADE algorithm proposed in this study.

Research
Robot arm

Conclusions
In this research, inverse kinematics problem for a 7 degree of freedom serial robot manipulator was implemented to prove the accuracy and efficiency of the self-adaptive control parameters in Differential Evolution (ISADE) and the ISADE algorithm with searching space improvement (Pro-ISADE) algorithm. To evaluate the effectiveness of the two algorithms above, the results obtained from the ISADE algorithm as well as Pro-ISADE were compared with the results from the PSO (Pro-PSO) and DE (Pro-DE) algorithm. Experiments were performed with three Scenarios. In the first Scenario, an endpoint in the workspace is randomly selected. The purpose of this Scenario is to compare the convergence speed of the three algorithms. In the second Scenario, algorithm was used to calculate inverse kinematics of the robot for 100 points randomly selected in the working space. The aim of this Scenario 2 is to test the accuracy and efficiency of the algorithm when the end effector started at the same position, it went to any point in working space. Meanwhile, in the third Scenario, the algorithms solved the inverse kinematics problem when the end effector of the robot moved point to point that are located on a spiral trajectory in the workspace. The implementation experiments have shown, the ISADE algorithm gave much better results than other algorithms in term of: accuracy, execution time and number of generation. Besides, by improving the searching boundary for joints' variable, the Pro-ISADE, Pro-DE and Pro-PSO also improve the accuracy as well as processing speed and especially the quality of the value of the joints variable compared to the ISADE, DE and PSO, respectively. These optimal joints' values ensure the feasibility of the dynamic and control problem in the future. In short, with ISADE algorithm as well as Pro-ISADE, they have handled the inverse kinematic requirement very effectively both in term of accuracy and computation time. The Pro-ISADE algorithm not only improves the above two factors, but also improves the quality of the joints' variables.

Conflict of interest
"The authors declare no conflict of interest."