Dimensional Synthesis for Multi-Linkage Robots Based on a Niched Pareto Genetic Algorithm

The dimensional synthesis of multi-linkage robots has great significance for improving flexibility and efficiency. With the increase of the degree of freedom and restrictions on special occasions, the solution of dimensional synthesis becomes complicated and time-consuming. Theory of workspace density function, maneuverability, and energy expenditure had been studied. With high flexibility and low energy consumption as the design goal, the method for dimensional and joint angle synthesis of multi-linkage robots was proposed based on a niched Pareto genetic algorithm. The Pareto solution set has been obtained. The method was verified by two application examples, which is occlusion of the solar salt evaporation pool and the secondary scattering of solid 2,2′-azobis(2,4-dimethylvaleronitrile). Through the application of NPGA (niched Pareto genetic algorithm) compared with KPCA (kernel principal component analysis), it can save 12.37% time in occlusion of one evaporating pool and reduce energy consumption by 3.85%; it can save 9.96% time in scattering of remain materials per barrel and reduce energy consumption by 1.77%. The study reduces the labor intensity of manual workers in the salt making industry, ensures the safe production of dangerous chemicals, and provides new ideas and methods for the dimensional synthesis of multi-linkage robots.


Introduction
A multi-linkage robot is a linkage mechanism composed of multiple movable joints connecting rigid rods in series [1,2]. It has the advantages of large range of working space, high flexibility, simple mechanical structure, good obstacle avoidance performance, and better avoidance of singularity of mechanism. They are widely used in medical surgery, underwater detection, construction, and many other engineering fields.
The fact shows that dimensional synthesis for multi-linkage robots has been paid much attention by scholars. Gosselin proposed a global condition index (GCI) that combines the working space and the Jacobian condition number to design the structure size of multi-linkage robots [3]. Based on the optimization goal of the subjective minimum robot criterion, Amir optimized the robot size with the Lagrange multiplier method [4]. José-Alfredo et al. used the means of a multi-objective optimization to carry out dimensional synthesis of a spherical parallel manipulator, the method of multi-objective evolutionary algorithm based on decomposition (MOEA/D) was proposed, objective functions were the dexterity, the variability of the dexterity, the maximum required torque, and the stiffness. The optimization of robots' workspace was performed using non-dominated sorting genetic algorithm-II (NSGA-II), multi-objective particle swarm optimization algorithm (OMOPSO), MOEA/D separately, the results showed better performance of MOEA/D algorithm [5]. Ramon et al. obtained the best linkage by means of hybridizing a local search (LS) approach and an evolutionary algorithm (EA).
The rest of the paper is organized as follows. In Section 2, three objective functions are researched by calculating appropriate cases, then the relevant calculation formula are derived. In Section 3, the parameters and design process of NPGA are determined. In Section 4, the mathematical model of operation is established. In Section 5, two application examples of parameter synthesis are given to verify the proposed method. Finally, we present our conclusions in Section 6.

Objective Functions of NPGA
NPGA adopts Pareto winning relation to select tournament. We choose objective functions combining the application of the n-R linkage robot based on the fitness sharing mechanism.

Workspace Density Function
The mechanical mechanism workspace is divided into equal-sized units uniformly. The ratio of point number in each unit to point number in entire workspace end-effector can reached divided by unit volume is workspace density [29]. The greater workspace density, the better flexibility of the mechanism. The Monte Carlo method conducts statistics and tests with random numbers and takes the mean value and probability as the numerical solutions to unsolved problems, also known as Monte Carlo simulation. We study the workspace density, use the Monte Carlo method and then describe the shape, characteristics, and properties of flexible workspace more intuitively with the help of computer technology in the view of probability statistics. In the paper, calculation and simulation are implemented using the method combined with density function.

Calculation of Monte Carlo method
The principle of the Monte Carlo method for workspace solution is where W is workspace, k is the generalized joint variable, M(k) is forward kinematics expression, and G is the joint space sets. The Monte Carlo method generates a large number joint space points randomly, map them to flexible workspace, obtain corresponding end-position data, and then computes the boundary of the flexible workspace statistically. The specific flow of the algorithm is as follows: 1.
Divide the samples of each joint angle, construct the sample space G according to the definition of flexible workspace; 2.
Calculate the positions end-effector reached on the basis of forward kinematics expression M(k); 3.
Construct the flexible workspace W using outputted the end-effector position. The larger the sample size, the geometric contour will be clearer and the flexible workspace will be more accurate.
The forward kinematics equation was established based on Denavit-Hartenberg method (D-H method), according to the geometric relationship between the robot end-effector position, bar lengths, and joint angles in Figure 1.
Therefore, the P position vector of the end-effector can be expressed by Equation (6).
The forward kinematics equation in rectangular coordinate system and in polar coordinate system can be expressed in Equations (9) and (10) respectively.

Construction of Density Functions
Monte Carlo method calculates the workspace of series manipulator. According to the characteristics that the workspace density function conforms to the random density function, the mathematical model is built to describe the workspace density.
Assume that f i (p) is the density function of the manipulator bars, which p = (R(x), t) ∈ G. The workspace density function of nR linkage manipulator can be expressed as: Among them, * represents the convolution operation on G; represents the convolution operation on R n .
The random variables can be continuously evaluated in the workspace which is a continuous space. Therefore, all the continuous random variables can be expressed by the probability density function, which is shorted PDF. The multivariate normal distribution (Gaussian distribution) density function is defined as Equation (13) considering that x is a multidimensional variable.
Strictly generalization of Equation (13) obtain Equation (14), where µ is mean value, is order symmetric positive definite matrix.
The solution of µ, are as follows in Equations (15) and (16) respectively.
Since the sum of discrete probability distributions is 1, integration of probability density function is 1 always, as shown in Equation (17).
The continuous probability density function is not limited to 1, so in the study of workspace, all continuous random variables are measurable, and all continuous distributions have density. For a NR manipulator, if each link angle has M states, the samples number is K = M n in the linkage manipulator workspace, the variance of this sample space can be written as Equation (18).
When M is sufficiently large, the variance of the actual workspace density function can be rewritten as Equation (19) lim where 1,2,...,n = n i=1 σ 1 2 A n , A n is n order identity matrix, σ 1 2 is variance. Therefore, according to convolution iteration and extreme value theorem, the workspace density function constructed by Gaussian probability-density function is expressed in Equation (20), µ is expected value.
All angles have no limit, 20,000 sample points are taken to calculate the workspace nephogram based on Monte Carlo method, the results are shown in Figure 4.        All angles have no limit, we also take 20,000 sample points, the workspace nephogram based on Monte Carlo method are shown in Figure 6. The four views comprehensively describe the workspace density by the Monte Carlo method.
Algorithms 2020, 13, x FOR PEER REVIEW 9 of 22 All angles have no limit, we also take 20,000 sample points, the workspace nephogram based on Monte Carlo method are shown in Figure 6. The four views comprehensively describe the workspace density by the Monte Carlo method.

Maneuverability
When determining the parameters of the linkage robot, flexibility is an important considered factor. Yoshikawa [30] defined the maneuverability and gave the solution method. The flexibility variation of the linkage robot is demonstrated accurately by quantitative value. The maneuverability

Maneuverability
When determining the parameters of the linkage robot, flexibility is an important considered factor. Yoshikawa [30] defined the maneuverability and gave the solution method. The flexibility variation of the linkage robot is demonstrated accurately by quantitative value. The maneuverability can be solved by Equation (21).
Among them, Ψ describes the workspace flexibility of the linkage robot, J is Jacobian matrix. Therefore, the Jacobian matrix plays an important role in maneuverability solution. Combined with Equation (6), the relationship that between the end-effector position X, Y and the angles θ 1 , θ 2 , θ 3 in Figure 1 is shown in Equation (22).
Take the partial derivative of that, Equation (23) is gained.
We define the Jacobian matrix, as shown in Equation (24), each term can be gotten from Equation (25).
The maneuverability of the n-R linkage robot can be written as Equation (27).

Energy Expenditure
In Figure 1, define the quality of l 1 , l 2 , l 3 is m 1 ,m 2 ,m 3 respectively. Without considering linkage robot's friction, energy expenditure includes kinetic energy and potential energy. In this study, energy expenditure is analyzed taking 3R linkage robot as an example to obtain the energy expenditure of n-R linkage robot. The kinetic energy of l 1 is shown in Equation (28). The kinetic energy of l 2 is shown in Equation (29). The kinetic energy of l 3 is shown in Equation (30).
The total potential energy of l 1 , l 2 , l 3 , is represented in Equation (31). So far, the energy expenditure of the 3R linkage robot can be expressed in Equation (32).
The following is energy expenditure calculation of the n-R linkage robot. The centroid position of the rod is provided in Equation (33).
The centroid velocity can be expressed with Equation (34).
The kinetic energy of l n can be derived as Equation (35).
The potential energy of l n can also be written as Equation (36).

Niche Pareto Genetic Algorithms
Based on Darwinian evolution theory, evolution is in the process of natural selection and adaptive. For intractable nonlinear complex problems in traditional algorithms, it has important application value that design and construct new evolutionary algorithm to solve the multi-objective optimization problem. We gain one solution in Pareto sets usually in traditional algorithms. However, more-or even all-Pareto sets can be solved using the genetic algorithm.
The pennant selection system of NPGA is different from comparisons among individuals distinctly. NPGA also selects other individuals for comparison. Let individuals evolve in certain living environment for forming multiple niches, then individuals will always evolve in each niche until they reach a peak in the niche. Therefore, the global optimal solutions are obtained [31].

NPGA Sharing Function
In order to ensure data diversity, according to individual similarity, population is divided into different sub-populations in parameter space and gene space. The sharing function between populations can express in Equation (39) where d ij is the distance between population i and j, σ share is a shared parameter controlled allowable sharing degree, α = 1 generally.

NPGA Fitness
The fitness between populations is adjusted with sharing function of populations' similarity degree. In that way, the sub-population is divided by new fitness in iterative evolution, so the niche evolution environment [32] is produced. The population fitness can be expressed as Equation (40) f where M is the adjacent population number of the ith population, f si is fitness after sharing, f i is fitness before sharing, and p is population size.

Choice of Coding Method
Each gene value is represented by floating point numbers in the algorithm, the encoding length is the number of variables. The method of floating-coding can deal with complex variable constraint conditions. It is convenient for genetic searching. Therefore, the floating-coding method is selected for the established mathematical model the variables of which are length of a structural member and joint angle. The mathematical model's aims are the density function, maneuverability, and energy expenditure. At the same time, it is also easy to mix NPGA with classical computing methods. In the floating-point encoding, The symbol bit is represented by S: 0 is positive, 1 is negative; N represents the binary bits of an int; The exponent bit is denoted by E; Mantissa is represented by M;

Choice of Operators
The purpose of choosing appropriate operators is to ensure the computational efficiency and improve the convergence of the algorithm. The survival of the fittest is carried out with elitist strategy, the worst individual is always replaced by the best, and it can limit the reproduction of certain individuals into the next generation. In this paper, the niche technique and Pareto optimal individual are applied to solve problems such as calculating the optimal solution of individual sub-objective function. The number of niches (n i ) is the number i of individuals in the niche segment of the genotype space (R n ). If two non-inferior individuals a, b are selected based on the niche, let n i (a) and n i (b) are respectively the niche number of individuals a, b in the niche i. Then, the corresponding individual is selected from mini(n i (a), n i (b)) (the individual is selected with the smallest number of niche).

Choice of Crossover Operator
Crossover operation is the main method to generate new individuals, which determines the global search ability of the algorithm. It should not damage the good character of individual coding string too much, but also effectively generate good new individuals. Adopting the method of uniform block crossover operator, the operation between individual X 1 and X 2 create two new individuals, as shown in Equation (41), where τ is a constant, and such a crossover operation is called uniform arithmetic crossover.

Choice of Mutation Operator
The mutation operation determines the local searching ability of genetic algorithm. Its purpose is to change the local searching ability of genetic algorithm and maintain the diversity of population. Because of the choice of floating-point encoding, a Gaussian mutation operator is used in this study. We define the value range of aberrance point that is x n ∈ [Q n min , Q n max ], then mean value µ and variance σ 2 with normal distribution can be expressed as Equations (42) and (43) separately.
Then new genic value may be determined in Equation (44), in which R is the sum of Uniformly distributed random Numbers in [0,1].

Operation Flow Chart
The operation flow chart of NPGA is shown in Figure 7. The procedure is as follows: 1.
Initialization of algorithm, produce new initial population and set up genetic parameters; 2.
Calculate the sharing degree;

Calculate individuals' fitness of population according to the individual sharing; 4.
Process of selecting, crossing, and mutating; 5.
Comparing of fitness size between offspring and last generation individuals; 6.
Offspring replace last generation individuals, forming a new population; 7.
The algorithm stops until the termination condition is activated, otherwise go back to Step 2.

Mathematical Model of Dimensional Synthesis
There are many performance indexes of a multi-linkage robot, especially on specific occasions and in specific engineering applications. Workspace density function, maneuverability, and energy expenditure are ordinary. Determining the parameters of the linkage robot to make it have higher flexibility and lower energy consumption has important application significance. In this study, we define the parameter variable is Equation (45).

Mathematical Model of Dimensional Synthesis
There are many performance indexes of a multi-linkage robot, especially on specific occasions and in specific engineering applications. Workspace density function, maneuverability, and energy expenditure are ordinary. Determining the parameters of the linkage robot to make it have higher flexibility and lower energy consumption has important application significance. In this study, we define the parameter variable is Equation (45).

Engineering Example of the Method Application and Analysis
The method of dimensional synthesis can be verified better in actual engineering. In this paper, the method is applied in salt industry and chemical industry. In addition, kernel principal component analysis (KPCA) was used to calculate the two-link and three-link robots in the study and application examples [33], and the results were compared with those obtained by NPGA in this paper.

Method Application in Salt Industry
Solar salt is the process of making salt from seawater, including moisture absorption, halogen generation, crystallization, and other processes. The salt-making process of seawater is generally stored in the evaporation pond for 7-15 days. According to the influence from weather, brine concentration, and other external factors, occlusion of the solar salt evaporation pool is necessary, it has great significance for improving the solar salt efficiency. Without manual retracting and retracting, the shielding cloth (as shown in Figure 8) will get off the track and become entangled, affecting the evaporation of the salt field. At present, the occlusion cloth is put in with the help of human, as shown in Figure 8 Due to the large number of evaporation pools, so it is very labor-intensive, time-consuming, and also less efficient.

Engineering Example of the Method Application and Analysis
The method of dimensional synthesis can be verified better in actual engineering. In this paper, the method is applied in salt industry and chemical industry. In addition, kernel principal component analysis (KPCA) was used to calculate the two-link and three-link robots in the study and application examples [33], and the results were compared with those obtained by NPGA in this paper.

Method Application in Salt Industry
Solar salt is the process of making salt from seawater, including moisture absorption, halogen generation, crystallization, and other processes. The salt-making process of seawater is generally stored in the evaporation pond for 7-15 days. According to the influence from weather, brine concentration, and other external factors, occlusion of the solar salt evaporation pool is necessary, it has great significance for improving the solar salt efficiency. Without manual retracting and retracting, the shielding cloth (as shown in Figure 8) will get off the track and become entangled, affecting the evaporation of the salt field. At present, the occlusion cloth is put in with the help of human, as shown in Figure 8 Due to the large number of evaporation pools, so it is very laborintensive, time-consuming, and also less efficient. As shown in Figure 9, a 2R linkage robot is designed instead of manually retracting or releasing of the shielding cloth. The two-lever robot actuator can clamp the blocking cloth to ensure the smooth rolling and releasing of the blocking cloth in the salt field, so as to save labor and improve the rolling and releasing efficiency of the blocking cloth in multiple evaporation pools in the salt field. As shown in Figure 9, a 2R linkage robot is designed instead of manually retracting or releasing of the shielding cloth. The two-lever robot actuator can clamp the blocking cloth to ensure the smooth rolling and releasing of the blocking cloth in the salt field, so as to save labor and improve the rolling and releasing efficiency of the blocking cloth in multiple evaporation pools in the salt field. Algorithms 2020, 13, x FOR PEER REVIEW 16 of 22 Figure 9. Robot assistance.
The dimensional synthesis data of two-link robot obtained by the KPCA is shown in Table 1. after radian convert to angles. The design solution is based on NPGA method. For the calculation model constructed in this paper, the population size is set as 100 and the chromosome length is set as 32 bits in binary order to obtain the global optimal solution quickly based on experience. The choice of crossover operator affects whether it can converge to the optimal solution region, so the crossover operator is 0.08. If the mutation operator is selected larger, on the one hand, increasing the diversity of samples can also cause instability. On the other hand, if the selection is small, it can be difficult to find the global optimal solution. Thus, mutation operator is set to select 0.65, termination of algebra for 100 generations, and the result is shown in Figure 10. According to the data, the convergence is good, the global Pareto solution is set, some data are shown in Table 2. The dimensional synthesis data of two-link robot obtained by the KPCA is shown in Table 1. In Table 1, dimensional parameters, we take l 1 = 12.0000, l 2 = 16.0000, θ 1 ∈ [0, 36 • ], θ 2 ∈ [0, 120 • ] after radian convert to angles.
The design solution is based on NPGA method. For the calculation model constructed in this paper, the population size is set as 100 and the chromosome length is set as 32 bits in binary order to obtain the global optimal solution quickly based on experience. The choice of crossover operator affects whether it can converge to the optimal solution region, so the crossover operator is 0.08. If the mutation operator is selected larger, on the one hand, increasing the diversity of samples can also cause instability. On the other hand, if the selection is small, it can be difficult to find the global optimal solution. Thus, mutation operator is set to select 0.65, termination of algebra for 100 generations, and the result is shown in Figure 10. According to the data, the convergence is good, the global Pareto solution is set, some data are shown in Table 2.
It can be seen from the data comparison in Table 3 that the multi-linkage robot is used to assist the rolling and releasing of the shielding cloth. From the rolling and releasing process of 3-6 times per day, it takes time 6 12 for each evaporation pool under the KPCA method. While for the NPGA method, it takes time 5 26 . Only one evaporation pool saves 12.37% the rolling and releasing time, which improves the rolling and releasing efficiency of the evaporation pool shielding cloth in the salt making process. The energy consumption of multi-linkage robot was estimated from the change of joint angle, and the NPGA method reduced the energy consumption by about 3.85% compared with the KPCA method.  It can be seen from the data comparison in Table 3 that the multi-linkage robot is used to assist the rolling and releasing of the shielding cloth. From the rolling and releasing process of 3-6 times per day, it takes time 6'12'' for each evaporation pool under the KPCA method. While for the NPGA method, it takes time 5'26'' . Only one evaporation pool saves 12.37% the rolling and releasing time, which improves the rolling and releasing efficiency of the evaporation pool shielding cloth in the salt making process. The energy consumption of multi-linkage robot was estimated from the change of joint angle, and the NPGA method reduced the energy consumption by about 3.85% compared with the KPCA method.

Method Application in Chemical Industry
2,2′-azobis(2,4-dimethylvaleronitrile)(shorted ADVN) (molecular formula: C14H24N4) is widely used in industrial production for polyvinyl chloride(PVC), polyacrylonitrile (PAN), polymethyl methacrylate(PMMA), anticorrosive paint, etc. Due to its strong toxicity, flammability, and explosive risk, the workers should have little or no contact with it as much as possible. In this, the use of linkage robots is very important to ensure human health. In the process of barreled solid

Method Application in Chemical Industry
2,2 -azobis(2,4-dimethylvaleronitrile)(shorted ADVN) (molecular formula: C14H24N4) is widely used in industrial production for polyvinyl chloride(PVC), polyacrylonitrile (PAN), polymethyl methacrylate(PMMA), anticorrosive paint, etc. Due to its strong toxicity, flammability, and explosive risk, the workers should have little or no contact with it as much as possible. In this, the use of linkage robots is very important to ensure human health. In the process of barreled solid ADVN scattering, there is some surplus in the barrel. In order not to damage the material barrel, we designed the 3R linkage robot for the secondary scattering of remaining materials, as shown in Figure 11.
In this application, the height of the barrel is 59 cm, the maximum diameter is 40 cm, and bunghole diameter is 26 cm. The dimensional synthesis data of three-link robot obtained by the KPCA is shown in Table 4.
Design solution is carried out based on the NPGA method. According to the mathematical model constructed and the application example of two-link robot in salt field, population size is 150, chromosome length in binary is 32 bits, crossover operator is 0.07, mutation operator is 0.70, termination generation is 150. The output results are indicated in Figure 12. It is observed that data convergence is good. We gained the global Pareto sets, some data are listed in Table 5. ADVN scattering, there is some surplus in the barrel. In order not to damage the material barrel, we designed the 3R linkage robot for the secondary scattering of remaining materials, as shown in Figure  11. In this application, the height of the barrel is 59 cm , the maximum diameter is 40 cm , and bunghole diameter is 26 cm . The dimensional synthesis data of three-link robot obtained by the KPCA is shown in Table 4. after radian convert to angles. Design solution is carried out based on the NPGA method. According to the mathematical model constructed and the application example of two-link robot in salt field, population size is 150 , chromosome length in binary is 32 bits , crossover operator is 0.07 , mutation operator is 0.70 , termination generation is 150 . The output results are indicated in Figure 12. It is observed that data convergence is good. We gained the global Pareto sets, some data are listed in Table 5.  Links     Links In Table 5, synthesizing various parameters, we take l 1 = 16.0000, l 2 = 22.0000, l 3 = 26.0000, after radian convert to angles.
From the data comparison in Table 6, the average processing time of a single barrel of ADVN in solid state is time 4 11 spent under the KPCA method, while that under the NPGA method is time 3 46 saved by single barrel processing. It improves the processing efficiency of powder making and reduces the manual time spent. From the change of joint angle, the NPGA method can reduce the energy consumption about 9.96% compared with the KPCA method.  Figure 13 shows the remaining ADVN in the barrel. Figure 14 shows the scattered ADVN that it has been completely separated from the barrel. This meet the requirement with the next step into the grinder.

Conclusions
A component's size is an important issue in mechanical design. Reasonable selection of construction size is of great significance for reducing energy consumption and improving the performance of multi-linkage robots. Based on theory and experiment, this paper studies the calculation method of dimensional synthesis of multi-linkage robots and draws the following conclusions: • Based on the study of workspace, maneuverability, and energy expenditure, the NPGA method for dimensional synthesis of multi-linkage robot was proposed and applied. Then, the NPGA was applied; • The superiority of NPGA method was verified by comparing with the KPCA method in two applications. The study provided new ideas and methods for the salt field's unmanned mechanized production and design of hazardous chemical processing equipment; • For the method of dimensional synthesis of multi-linkage robots, it is possible to obtain a locally convergent solution set under a special environment and specific constraint conditions, but the Pareto solution set with global convergence cannot be obtained. It is possible to obtain a locally convergent solution set, and the algorithm needs to be further improved.

Conclusions
A component's size is an important issue in mechanical design. Reasonable selection of construction size is of great significance for reducing energy consumption and improving the performance of multi-linkage robots. Based on theory and experiment, this paper studies the calculation method of dimensional synthesis of multi-linkage robots and draws the following conclusions: • Based on the study of workspace, maneuverability, and energy expenditure, the NPGA method for dimensional synthesis of multi-linkage robot was proposed and applied. Then, the NPGA was applied; • The superiority of NPGA method was verified by comparing with the KPCA method in two applications. The study provided new ideas and methods for the salt field's unmanned mechanized production and design of hazardous chemical processing equipment; • For the method of dimensional synthesis of multi-linkage robots, it is possible to obtain a locally convergent solution set under a special environment and specific constraint conditions, but the Pareto solution set with global convergence cannot be obtained. It is possible to obtain a locally convergent solution set, and the algorithm needs to be further improved.