Multi-objective optimization design of parallel manipulators using a neural network and principal component analysis

. In this work, a multi-objective optimization design method is proposed based on principal component analysis (PCA) and a neural network to obtain a mechanism’s optimal comprehensive performance. First, multi-objective optimization mathematical modeling, including design parameters, objective functions, and constraint functions, is established. Second, the sample data are obtained through the design of the experiment (DOE) and are then standardized to eliminate the adverse effects of a non-uniform dimension of objective functions. Third, the ﬁrst k principal components are established for p performance indices ( k < p ) using the variance-based PCA method, and then the factor analysis method is employed to deﬁne its physical meaning. Fourth, the overall comprehensive performance evaluation index is established by objectively determining weight factors. Finally, the computational cost of the modeling is improved by combining the neural network and a particle swarm optimization (PSO) algorithm. Dimensional synthesis of a Sprint (3RPS) parallel manipulator (PM) is taken as a case study to implement the proposed method, and the optimization results are veriﬁed by a comprehensive performance comparison of robots before and after optimization.


Introduction
Parallel manipulators (PMs) have become a research hotspot because of their excellent dynamic performance, simple inverse kinematics, and high stiffness (Yang et al., 2022). Sprint (3PRS) (Chen et al., 2014), Tricept (3UPS-UP) (Farooq et al., 2021), and delta (Wang et al., 2014) robots have been successfully commercialized, and the 3RPS PM is used as a parallel module in the five-axis machining center for blade surface machining (Arabshahi and Novinzadeh, 2015). Comprehensive performance of PMs is crucial for machining precision. Therefore, multi-objective optimization is necessary to achieve optimal comprehensive performance.
PMs have smaller workspaces compared with serial robots. Therefore, maximizing the workspace volume, especially the regular workspace volume, is a major task of a Since multiple performance indices are often presented as competitive coupling relations, this makes it difficult to obtain an optimal solution with a definite physical meaning. The comprehensive objective method is most commonly used to fix this problem, and it defines the weighted sum of objective functions as a comprehensive index Zhang and Nelson, 2011). Xu et al. (2017) optimized the hybrid mechanism's comprehensive performance through a comprehensive objective method wherein the workspace volume and stiffness performance were considered as objective functions. Huang et al. (2018) adopted this method to optimize reconfigurable a PM's workspace volume, global dexterity, and stiffness indices. Various combinations of weighting factors are provided to decisionmakers. However, the weight coefficients in this method are often determined artificially, and the dimension and unit of the objective function are inconsistent, which leads to the limitation of the application of this method (Fan et al., 2022). Yang et al. (2019) derived a globally optimal comprehensive performance by proposing a multi-objective optimization game algorithm. Globally optimal solutions are derived based on their weight factors rather than an objective function weight sum. The method is simple but needs to manually set the weight factors and calculate the objective function's values at each discrete node in advance. Discrete node density decides the precision. As the density of discrete nodes increases, the calculation cost increases exponentially. For instance, for a PM with 3 degrees of freedom (DOFs), a population n p , and generations n g , if the node density in the workspace is increased to 2 times that of the original, the computational cost will be n p n g 2 3 times that of the original (Yang et al., 2019;Wu, 2017).
The Pareto front method has become a research hotspot of multi-objective optimal design methods because it can obtain multiple sets of non-dominated solutions (Altuzarra et al., 2011;Ling et al., 2022). Cui et al. (2015) combined a genetic algorithm (GA) and insight software to optimize the kinematics and kinetic flexibility of the tunnel tube assembly system's 3-DOF PMs and obtained multiple sets of nondominated solutions. Qi et al. (2018), Sun and Lian (2018), and Sun et al. (2019) defined the sum of minimum squares (i.e., the cooperative balance point) from the non-dominated solutions as the global optimal solution. However, the physical explanation for this method is not clearly defined.
Establishing workspaces' analytic formulas is challenging, especially for regular workspaces, so the GPIs often obtained in terms of numerical methods include the Monte Carlo method (Wang et al., 2015) and the discretenode method (Xu et al., 2018). High computational costs greatly reduce the optimized design's computational efficiency. The response surface model (RSM), including a back-propagation (BP) neural network (Zhang and Gao, 2008;Gao et al., 2010), multivariate regression (Yang et al., 2022;Lian et al., 2019), and Gaussian regression (Vinh et al., 2019), is often used to fix this problem, wherein the BP neural network has attracted extensive attention due to its strong learning and nonlinear fitting capabilities. Therefore, the BP neural network model is adopted to establish the RSM of the GPIs and the design parameters in this work.
The main contribution of this work is to propose a multiobjective optimization design method in terms of factor analysis, BP neural network, principal component analysis (PCA), and particle swarm optimization (PSO) algorithm (Sun and Lian, 2018; to achieve the optimal comprehensive performance of the mechanism. The proposed method has a clear physical explanation and a low computational expense.
The remainder of the work is arranged as follows. The structure description of the 3RPS PM is presented in Sect. 2. Objective functions of the 3RPS PM are presented in Sect. 3. The multi-objective optimization design of the 3RPS PM is introduced in Sect. 4, and the optimization results are discussed in detail. Finally, conclusions are drawn in Sect. 5.

Procedure of the multi-objective optimization design method
The basic idea of the multi-objective optimization design for the PMs proposed in this work is to first use PCA to construct an overall comprehensive performance evaluation index with objectively determined weight factors and clear physical significance. Then a response surface model between objective functions and design parameters is constructed based on a neural network to solve the problem of the high computational cost of the objective functions in the workspace. Figure 1 shows the procedure of the proposed multi-objective optimization design method of PMs based on a neural network and PCA. The steps are described as follows.
Step 1. Establish the multi-objective optimization mathematical model, including design parameters, objective func-tions, and constraint functions, wherein the limitations of the design parameters, actuation lengths, and joint angles should be given.
Step 2. Obtain the sample data based on the design of the experiment: Latin hypercube sampling (LHS) is recommended in this work, and then it is necessary to standardize the original sample data to eliminate adverse influences caused by non-uniform dimensions of objective functions. The Z-score method is adopted in this work.
Step 3. Establish the correlation matrix of multiple performance indices and the first k principal components according to the cumulative contribution rate (the ratio of the sum of the first k eigenvalues of the correlation matrix to the sum of all the eigenvalues). The i principal components can be established based on the corresponding eigenvectors of the correlation matrix. Finally, the overall comprehensive performance index is defined as the sum of the product of the first k eigenvalues and the corresponding principal components.
Step 4. Establish the RSM between GPIs and design parameters based on the neural network to avoid the high computational cost of the objective functions in the workspace.
Step 5. Combine the overall comprehensive performance index and RSMs to carry out the optimization design of PMs to achieve optimal design parameters with a clear physical explanation. Figure 2 shows a spatial 3RPS PM with 3 DOFs. Three identical branches connect the base and mobile platform at points B i and A i with revolute joints and spherical joints, respectively. An upper rod and a lower rod connected by a prismatic joint constitute limbs. Inputs to the mechanism refer to the three prismatic joints. Coordinate frames B i {x i , y i , z i }, o {x, y, z}, and O {X, Y, Z} are attached to the limb, moving platform, and base, respectively (see Fig. 2). Both A 1 A 2 A 3 and B 1 B 2 B 3 are equilateral triangles with |oA 1 | = |oA 2 | = |oA 3 | = r 2 and |OB 1 | = |OB 2 | = |OB 3 | = r 1 . Each of the revolute joint axes is perpendicular to vector OB i and lies on the X-Y plane. Geometric and physical parameters before optimization are given as moving platform height h = 0.05 m, r 1 = 1.7 m, elasticity modulus E = 200 GPa, r 2 = 0.8 m, shear modulus G = 77 GPa, and d i = 0.4 m (diameter of the link). The 3RPS PM's multi-objective optimization aims to select appropriate r 1 and r 2 values to achieve maximal comprehensive performance.

Objective functions of the 3RPS PM
The 3RPS PM is to be used as the parallel module of the high-speed machining equipment, wherein dynamic performance, motion or force transmission performance, and operation workspaces strongly influence the mechanism's comprehensive performance. Therefore, the regular workspace volume, kinematic performance, and natural frequency are considered the three objective functions in this paper.

Workspace evaluation
A PM's workspace is smaller than those of serial mechanisms. Thus, PM optimization should primarily maximize workspace volume. A regular workspace volume (Yang et al., 2019) is recommended as the workspace evaluation index due to some performance indices being poor at the boundary of the reachable workspace, and the reachable workspace does not apply to trajectory planning (Babu et al., 2016). The method proposed in Yang et al. (2019) is adopted to calculate the regular workspace volume, the reachable workspace is first classified into n layers with an equal interval height of d z , and then the inscribed circle of the workspace boundary in each layer is calculated. The 3RPS PM's regular workspace volume is expressed as follows: where A i = 2π ρ 2 i is the inscribed circle's area in the ith layer. ρ i is the radius. The minimum radius of inscribed circles is set to 50 mm in the work to guarantee the necessary operation workspace of the 3RPS PM.

Kinematic performance evaluation
The kinematic condition number is often adopted as the kinematic index of PMs to measure their dexterity, manipulability, and singularity. The main limitation of the condition number index is the non-uniform dimension of the Jacobian matrix for the PMs with mixed DOFs. Characteristic length is often adopted to fix this problem. However, due to the subjectivity of characteristic length selection and the fact that only some elements of the Jacobian matrix are involved in characteristic length calculation, the physical explanation for the characteristic length is ambiguous (Zang et al., 2019;Wang et al., 2019). The motion or force transmission index (Wu et al., 2023;Liu et al., 2012) is recommended in this paper due to its advantages of being frame-free and dimensionless.
The ith limb's twist system is derived through screw theory: where S ij is the j th twist of the ith limb, l i is the unit vector along the vector B i A i . R i represents a rotation matrix from limb coordinate frame {B i } to global coordinate frame {O}. e i in Eq. (2) is presented in e 1 = 1 0 0 T , e 2 = 0 1 0 T , and e 3 = 0 0 1 T .
The transmission wrench of the ith limb can be obtained based on its reciprocal product with the twist (S ij ) of the passive joints being zero: The ith limb's input twist screw is All of the inputs except the one in the ith limb are locked to derive the instantaneous 1-DOF output twist screw S Oi : where U i is the mechanism's constraint wrench system. Accordingly, the ith limb's transmission and input indices are expressed as follows: Both of the output and input indices are considered to obtain the local transmission index (LTI): The global transmission index (GTI) is obtained to determine the global kinematic performance in a PM's regular workspace: where N is the number of discrete grid nodes in the prescribed workspace.

Dynamic performance evaluation
The natural frequency is the frequency at which a system tends to oscillate in the absence of external forces or damping. The first natural frequency, i.e., the fundamental frequency (Ganesh and Rao, 2020), is recommended in this paper due to its comprehensive consideration of the stiffness, mass, and damping of PMs. The global independent generalized displacement coordinates (IGDCs)  are recommended as the elastodynamic modeling method of PMs in this paper due to the Lagrangian multipliers and kinematic constraint equations not being required.
The ith link's elastodynamic equation expressed in the limb coordinate frames is given as follows considering that the link is discretized into one element (the number of elements has little influence on the calculation accuracy of the fundamental frequency; Yang et al., 2021): T , and L f Bi ( L m Bi ) and L f Ai ( L m Ai ) are the forces (moments) exerted on points B i and A i , respectively. K e is the element stiffness matrix. M e is the element mass matrix. L in the upper-left corner indicates that the vector is expressed in the local frame system.
Equation (9) can be further presented in the global coordinate system: where With a rigid assumption of the moving platform, its elastodynamic equation in the global coordinate system is where  (Yang et al., 2019).
Revolute joints are used to connect the base to link B i A i based on MPC theory. The boundary conditions in the local coordinate frame are as follows: where L ϕ Bx i and L ϕ Bz i are point B i 's angular displacement coordinates along the x i and z i axes, respectively. Boundary conditions in Eq. (13) can be further expressed based on Ganesh and Rao (2020).
Accordingly, the global IGDC of the manipulator can be extracted as follows.
Finally, the elastodynamic control equation of the 3RPS PM is given by trix from the generalized displacement coordinate U i to the global IGDC U . It can be seen that the Lagrangian multipliers or the simultaneous kinematic constraint equations are not required in Eq. (16): all the kinematic constraints of the mechanism were considered in the global IGDC and the mapping matrix N i , which is beneficial for solving dynamic equations.
The angular frequency in radian per second is derived by The fundamental frequency is The PM's global dynamic index is

Multi-objective optimization analysis using the PCA and the neural network
The multi-objective optimization mathematical modeling of the 3RPS PM (Babu et al., 2013) is 0.5 ≤ r 1 ≤ 2.5, 0.5 ≤ r 2 ≤ 2.5, where L min = 1 m and L max = 4 m. L i is link B i A i 's length. ϕ i and θ i represent angles among the moving platform, base, and link, respectively, and θ max = ϕ max = 50 • . The PCA recombines the m dependent original indices into a new set of p (p < m) independent principal components (F i ) and determines an overall comprehensive index with a deterministic weight factor based on the variancebased method, which satisfies cov(F i , F j ) = 0. The physical meanings of the principal components are then explained by using the factor analysis method.
Consider the non-uniform dimensions of objective functions: it is necessary to convert original sample data X(X = [X 1 , X 2 , . . ., X p ], X i = [X 1i , X 2i , . . ., X ni ] T ) to the standardized X * before calculating the principal components. The Zscore method is adopted here.
where X * ij = X ij −µ j σ j , X * i = X * 1i , X * 2i , . . ., X * ni T , and σ j = var(X j ), µ j = E(X j ). Population characteristics are often represented by sample data in terms of the LHS method (Lian et al., 2017). Table 1 shows 276 sets of sample data designed based on the LHS method. Objective functions' standardized sample data of the 3RPS PM are presented through Eq. (21) (see Table 2). The correlation matrix R is established to reflect correlations of objective functions: where The correlation matrix can be defined in Eq. (23) for standardized data X * as follows: The correlation matrix of the objective functions obtained by Eq. (23) is shown in Table 3. Dynamic performance and motion or force transmission performance are negatively correlated with the workspace volume. Dynamic performance and kinematic performance present positive correlations, which is consistent with the results in Chen and Yang (2021). The first k principal components' cumulative contribution rate and the ith principal component are defined by where λ i and a i = (a 1i , a 2i , . . ., and a pi ) T are the ith eigenvalue and eigenvector of the correlation matrix, respectively. δ is the first k principal components' cumulative contribution rate. The first two eigenvalues' cumulative contribution rate is 92.65 %, which is greater than 85 % (Kuroda et al., 2011) (see Table 4). Objection function information is included in the first two principal components (see Eq. 25).
The physical explanation of principal components is given by introducing the factor analysis method. The factor-loading matrix is first given by The factor-loading matrix should be further simplified to explain principal components until each objective function has a large loading in only one principal component. The maximum orthogonal rotation method of square differences is usually adopted (Mardian et al., 1979): where is the orthogonal matrix. Repeat Eq. (28) until the allowable error is larger than the change in B's relative variances.
The kinematic performance and the dynamic performance have larger factor loadings for the first principal component (see Table 5). Therefore, the kinematic and dynamic performance factor is named for the first principal component and presents a PM's dynamic performance and kinematic performance. The workspace volume has a larger factor loading on the second principal component; thus, the workspace factor is named for the second principal component, which reflects the mechanism's workspace volume. So far, the physical explanation of the principal components has been fully defined through the factor analysis method.
Accordingly, the 3RPS PM's overall comprehensive performance evaluation index is Generally, the GPI is obtained by calculating the variance or mean of the function value of nodes in the regular workspace. Computational costs increase exponentially with increased discrete points. The high computational intensity reduces the optimization design efficiency. The RSM between design parameters and GPIs is established using the BP neural network in this work to improve the computational intensity for its powerful nonlinear fitting effect. The BP neural network's hidden layer sizes are set to 30. The recommended acceptance levels of RMSE and R 2 in Agarwal and Renaud (2004) are below 0.2 and above 0.9, respectively. Table 6 shows the optimal accuracy assessment for the BP-based mapping model through five tests that verify the accuracy and effectiveness of the BP-based RSM.
The population-and evolution-based PSO intelligence algorithm determines their search directions based on their ve- locities and replaces the GA's mutation and crossover. The PSO algorithm's convergence efficiency was proven in . Thus, the PSO algorithm is used in this work. Table 7 shows the PSO algorithm's setting that is recommended in Sun et al. (2018). The 59 iterations are performed to converge the objective function to 3.5945, and the optimal design parameters are r 1 = 1.4340 and r 2 = 1.0660 (see Fig. 3). The calculation time of the proposed algorithm is 192 s on a computer 3.00 GHz CPU, while the PSO algorithm is 6285 s and 96.94 % of the computational cost is saved. Optimized dynamic performance, kinematic performance, and workspace volume increase by 44.57 %, 6.74 %, and 283.69 %, respectively; after optimization, the overall comprehensive index F increases by 118.92 % (see Table 8).
The multi-objective optimization design results verified the effectiveness of the proposed method. Dynamic performance, kinematic performance, and regular workspace have been enhanced as shown in Table 8. The distributions of the performance indices are shown in Figs. 4-5, the maximum radius of the workspace has been increased from 0.096 to 0.156 m, the maximum LTI has been increased from 0.95 to 0.99, and the maximum first natural frequency has been increased from 28.42 to 66.57 Hz, which is consistent with the results in Table 8. It is worth noting that, since the mechanism has a relatively excellent kinematic performance before optimization, the kinematic performance has not been significantly improved after optimization.

Conclusions
The results of multi-objective optimization design significantly affect the comprehensive performance of PMs. Ob- taining the mechanism's optimal comprehensive performance with clear physical significance is challenging because of weight factors' subjectivity and multiple performance indices' non-uniform dimensions. Meanwhile, the high computational costs of the GPI undermine PMs' multiobjective optimization design efficiency. In this paper, the multi-objective optimization design method based on the PCA, factor analysis, PSO algorithm, and neural network was proposed to fix these problems, and the multi-objective optimization of the 3RPS PM was taken as a case study to implement the proposed method.
Considering the fundamental frequency, the motion or force transmission performance, and the regular workspace volume of the 3RPS PM to be objective functions, the results of PCA and factor analysis showed that the first principal component reflected the mechanism's dynamic and kinematic performances. The second one reflected the mechanism's regular workspace volume. The computational cost of the proposed approach was reduced by 96.95 % compared with the PSO algorithm. Meanwhile, the mechanism's overall comprehensive performance increased by 118.92 %, and the dynamic performance, kinematic performance, and workspace volume increased by 44.57 %, 6.74 %, and 283.69 %, respectively.
Although neural networks have a strong fitting ability, they have not revealed the inherent mechanism between the objective functions and the design parameters. In our future work, we will study the analytical formula between the objective function and the design parameters.    Code and data availability. The corresponding author can provide MATLAB and Ansys files to support the work. Review statement. This paper was edited by Zi Bin and reviewed by three anonymous referees.