# Multicriteria Optimization Design for End Effector Mounting Bracket of a High Speed and Heavy Load Palletizing Robot

Authors: Ying He; Jiangping Mei; Jiawei Zang; Shenglong Xie; Fan Zhang
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

End effector mounting bracket is an important load bearing part of high speed and heavy load palletizing robot, which is located at the most distant point in robot rotation radius and frequently works in complex conditions such as start-stop, switch direction, and acceleration and deceleration motion; therefore, optimizing design for its structure is beneficial to improve the dynamic performance of robotic system and reduce energy consumption. Firstly, finite element model of end effector mounting bracket was established, and its accuracy was verified by contrastive analysis of modal test result and finite element model. Secondly, through modal analysis, vibration response test, frequency response analysis, and the static analysis, taking inertia into account, the mass is minimized, the maximal stress is minimized, the maximal deformation is minimized, and the first natural frequency is maximized as the optimization objectives are determined; the design variables were selected by sensitivity analysis, taking their value range as the constraint conditions; approximation models of objective functions were established by the Box-Behnken design and the response surface methodology, and their reliability was validated; to determine weighting factor of each optimization objective, an analytic hierarchy process based on finite element analysis (FEA + AHP) method was put forward to improve the objectivity of comparison matrix; subsequently, the multicriteria optimization mathematical model was established by the methods mentioned above. Thirdly, the multicriteria optimization problem was solved by the NSGA-II algorithms and optimization results were obtained. Finally, the contrastive analysis results between optimized model and initial model showed that, in the case of the maximum stress and deformation within allowable values range, the mass reduction was 17.8%; meanwhile, the first natural frequency was increased, and vibration response characteristics of the entire structure were improved significantly. The validity of this optimization design method was verified.


---

## Body

## 1. Introduction

Because palletizing robots can accurately and efficiently take the place of humankind in handling changing palletizing tasks, they have been widely applied in logistics operations [1, 2]. With the accelerating pace of production and as product scale expands unceasingly, putting forward higher requirements for palletizing robot performance such as work efficiency, load capacity, and motion stability, this will unceasingly present new challenges to robot structural design and reliability. The end effector mounting bracket is an important part which is used to mount end effector and bear working load; it is located at the most distant point in robot rotation radius and frequently works in start-stop, switch direction, acceleration and deceleration motion, and other complex working conditions; therefore, optimizing the design for its structure and reducing its mass in the premise of guaranteeing the strength, rigidity, and vibration stability not only are beneficial to improve the dynamic performance of robotic system but also reduce energy consumption.Domestic and foreign scholars have done a lot of studies on optimal design of the robot structure. Kunpeng et al. researched structure optimization problem on ER300 type palletizing robot forearm and end bracket based on the finite element method (FEM), using minimal mass as the optimization objective, stiffness and strength as the constraints, and plates thickness and holes diameter as design variables [3]. Oral and Kemal Ider researched structure optimization problem on flexible robot arm based on the FEM. Using stress and deformation as constraints, section size of each robot arm as design variables, and the sequential quadratic programming algorithm as solving method optimized the structure of the robot arm [4]. Albers and Ottnad researched structure optimization problem on the arms of the humanoid robot ARMAR III based on the mechanical and electrical cosimulation optimization method, effectively reducing the mass of the arm [5]. Ghiorghe researched lightweight design on a RRR type industrial robot based on the FEM, used size parameters as design variables, and reduced the mass of robot [6]. Saravanan et al. researched lightweight design on 2-link and 3-link planar robots, using minimized torque and maximized manipulability measure of robot as objectives, maximized deformation and joint angle as constraints, and lengths and cross-sectional area parameters of mechanical arm as design variables [7]. Guo et al. researched lightweight design on manipulator elderly service robot based on the FEM, reducing the mass of manipulator under the condition of satisfying strength and stiffness requirements [8]. Yoshimura and Masui researched on integrated optimization of motion planning and structure of the robot, using the stiffness as constraint, by optimizing the section size of mechanical arm to reduce mass and inertia [9]. Anfu and Cheng optimized 6-DOF robotic welding big arm based on statics analysis via FEM; maximum stress and maximum deformation were decreased; the mass of the big arm was significantly reduced [10]. Rui and Ma optimized frame of aluminum ingot palletizing robot end effector based on FEM, under the condition of satisfying strength and stiffness requirements, reducing mass of the end effector frame [11]. Most robots involved in the above study work at lower speed and lighter load working condition; the corresponding research is mostly based on statics analysis; main factors considered are strength and stiffness; under the condition of satisfying strength and stiffness requirements, the mass of robot is reduced. However, with the constant improvement of the palletizing robot high speed and heavy load performance requirements, the influence of inertia force cannot be ignored, provided that the palletizing robot and the load are regarded as a multidegree vibration system; the system will be subject to a variety of exciting forces such as inertial force, joint torque, and heavy torque. In the working process of the robot, these exciting forces alter in real time, and their frequency increases with the increase of speed, but each order natural frequency of the palletizing robot and each posture do not alter with the increase of speed. When frequency of exciting force increases to reach natural frequency of the system, resonance may occur, affecting the working stability of palletizing robot, even causing damage of components. Therefore, for structure optimization design on high speed and heavy load palletizing robot, statics analysis is not enough; we should also do the modal analysis and use natural frequency of structure as constraint or optimization goal in the optimization design.In the other studies, carried out on statics analysis and modal analysis on the robot structure, Cao et al. [12] and Tian et al. [13], respectively, researched on the lightweight palletizing robot arm and reducer machine base of 7-DOF humanoid robot arm shoulder; modal analysis was simulated with FEM; however, the natural frequency was not considered in the optimization process.Liu et al. researched structure optimization problem on wafer handling robot arm, with the natural frequency of the robotic arm as the optimization objective, the end effector static deformation as constraints, and wall thickness and mass of the second and the third arm as optimization variables; the results show that optimization improved the natural frequency of the arm and reduced the amount of deformation of the end effector [14]; Ye et al. researched lightweight design on the humanoid robot frameworks. They adopted evolutionary structural optimization method and took into account some factors such as mass, maximum stress, and the first natural frequency [15]. However, in the above study, the basis of selection for the order of natural frequency was not discussed in detail, which has a certain degree of subjectivity and deficient guiding significance in practice, while the influence after the natural frequency improvement on the structure performance was not deeply discussed.In addition, the multicriteria optimization results are a Pareto optimal set, which need to set the weight coefficient of the optimization objectives to select the most satisfactory optimization design plan from the Pareto optimal set. Using analytic hierarchy process (AHP) obtains the weight of each performance index. AHP, which is a qualitative, quantitative, systematic, and hierarchical analytical decision-making method, was put forward by American operations researcher Saaty in the early 1970s [16, 17] and has been widely used to solve problems like performance evaluation, weight analysis, decision-making in education [18], management [19], energy [20], tourism [21], and other fields. However, applying AHP in the mechanical structure optimization study is relatively rare; in this paper, the AHP is introduced into the multicriteria structural optimization design. Additionally, the establishment of pairwise comparison matrix of the AHP usually relies on expert’s or designer’s experience or preferences, which has a certain degree of subjectivity [22–24]; it is difficult to obtain other researchers’ consensus. In this paper, an improved AHP based on results of finite element analysis is proposed in order to improve the objectivity of comparison matrix to a certain extent. The size of variable threshold of optimization objectives is obtained by the finite element analysis results. According to the size of variable threshold, the importance of optimization objectives can be obtained.In this paper, optimization object is end effector mounting bracket of MD-1200 type YJ palletizing robot. Firstly, the finite element model of end effector mounting bracket is established, and its veracity is verified by the contrastive analysis of modal test result and finite element calculating modal. Secondly, through modal analysis, vibration response test, frequency response analysis, and the static analysis, considering inertia, the mass is minimized, the maximal stress is minimized, the maximal deformation is minimized, and the first natural frequency is maximized as the optimization objectives are determined. Approximation models of optimization objective functions are established by the design of experiment (DOE) and the response surface method (RSM), and reliability of approximation models satisfies accuracy requirements through precision verification. We use the analytic hierarchy process based on finite element analysis (FEA + AHP) method to determine weighting factor of each optimization objective and establish multicriteria optimization mathematical model. Using the NSGA-II algorithm solves the multicriteria optimization problem and obtains optimization model parameters. Finally, it proves that the optimization design method is effective by contrastive analysis results of static analysis, model analysis, and frequency response analysis on the optimized model and initial model.


## 2. Multicriteria Optimization Design for the End Effector Mounting Bracket

The optimization object is the end effector mounting bracket of MD-1200YJ type palletizing robot. The MD-1200YJ type palletizing robot has 4 degrees of freedom; its 3D model is shown in Figure1. It is mainly composed of machine base, waist joint, small arm joint, drive arm of small arm, balancing weight, small arm drive link, attitude hold link 1, attitude hold tripod, attitude hold link 2, small arm link, big arm joint, big arm link, waist mounting bracket, end effector mounting bracket, and so on. The action of two parallelogram mechanical linkages (one consists of attitude hold link 1, big arm link, attitude hold tripod, and waist mounting bracket; the other consists of attitude hold tripod, attitude hold link 2, small arm link, and end effector mounting bracket) makes the undersurface of end effector keep parallel to the ground. The load capacity is 120 kg, maximum turning radius is 2400 mm, repeat positioning accuracy is ±0.4 mm, and waist maximum rotation speed is 85°/s, belonging to high speed and heavy load palletizing robot [25]. The end effector is mainly composed of AC servo motor, RV reducer, motor flange, end effector mounting bracket, some fasteners, and so on, with total mass of about 22.53 kg and mass of the end effector mounting bracket of about 9.79 kg, accounting for 43.5% of the total mass.Figure 1
3D model of MD-1200 type YJ palletizing robot.


### 2.1. Establishment of Finite Element Model of the End Effector Mounting Bracket

Before the finite element model is established, the 3D model is simplified by removing chamfer, fillet, screw holes, and others of the part. The simplified model is shown in Figure2. The material property and mesh setting are shown in Table 1 [26] and Table 2, respectively. The finite element model is shown in Figure 3.Table 1
Material property.
Item Property Material ZL110 Elastic modulus (N/m2) 6.89 × 1010 Poisson ratio 0.35 Shear modulus (N/m2) 3.19 × 108 Density (kg/m3) 2770 Tensile strength (N/m2) 1.65 × 108Table 2
Mesh setting.
Item Property Mesh type Tetrahedral Minimum element size (mm) 4.96532 Maximum element size (mm) 14.8961 Total number of nodes 19995 Total number of units 10884Figure 2
Simplified model of end effector mounting bracket.Figure 3
Finite element model of end effector mounting bracket.The following is a justification of the above mesh setting. First of all, element selection: tetrahedral mesh has been gaining popularity for solving finite element analysis with its simple, flexible, and strong adaptability to complex boundary (Chunge et al. [27]). The end effector mounting bracket is an irregular and complicated part; therefore, tetrahedral mesh is selected in this work. Second, number of elements: in order to ensure the accuracy of finite element analysis, the highest quality elements are adopted in the finite element software; therefore, the number of elements, total number of nodes, minimum element size, and maximum element size are obtained automatically.


### 2.2. Finite Element Model Accuracy Verification Based on the Modal Test

The modal test method is used to verify the accuracy of the finite element model. Experimental equipment includes LMS TEST.LAB vibration and noise analysis system, LMS SC305 type vibration and noise data acquisition system, three Lance LC0152T type unidirectional piezoelectric acceleration sensors, B&K Impact Hammer Type 8207, signal cable, IBM computer, crane, and elastic rope. We use single-point excitation and multipoint response test method that selects 21 pick points and 1 knock point (O point) on the part and then use impact hammer knock on the O point, respectively, along -x, +y, and +z directions; meanwhile, x, y, and z direction responses of every pick point (a total of 63 sets of data) are measured with three acceleration sensors, and measurement signals are acquired by the LMS SC305 system and are processed by the LMS TEST.LAB system; the test modal of the end effector mounting bracket is obtained afterwards. To measure the free modal of the part, it is overhung by a greatly flexible elastic rope to make it close to the free boundary conditions. The experimental environment is shown in Figure 4.Figure 4
The experimental environment.The free modal of the end effector mounting bracket is calculated by the finite element method. Comparing the FEM model with test modal for the first three modal shapes and the first five modal frequencies is shown in Figure5 and Table 3, respectively, where Ftestmodal is test modal frequency and Ffemmodal is finite element modal frequency.Table 3
Comparison between the test model and the FEM model for the first five modal frequencies (Hz).
Modal order F t e s t   m o d a l F f e m   m o d a l Relative error% 1 229.11 228.34 0.3404 2 448.41 449.75 0.3011 3 757.03 732.32 3.2628 4 789.54 768.58 2.6547 5 932.24 942.74 1.1263Figure 5
Comparison between the test model and FEM model for the first three modal shapes.
(a)
Comparison between the test model and FEM model for the first modal shape (b)
Comparison between the test model and FEM model for the second modal shape (c)
Comparison between the test model and FEM model for the third modal shapeThe comparison results indicate that the first three modal shapes of the calculating modal and test modal are consistent, and the relative error between the first five modal frequencies of the calculating modal and test modal is less than 5%. This shows that the accuracy of the established finite element model meets the requirements, and it can be used for the subsequent calculating simulation [28, 29].


### 2.3. Constraint Modal Analysis of the End Effector Mounting Bracket

According to the assembly relation of end effector mounting bracket and other parts, cylindrical surface constraints are applied in three bearing holes of its finite element model, and the first four modals are calculated by making use of Block Lanczos method. The first four natural frequencies are shown in Table4.Table 4
The first four natural frequencies of the end effector mounting bracket, Hz.
Order number 1 2 3 4 Frequency 409.57 802.07 1021.2 1312.1


### 2.4. Vibration Response Test

Using the vibration response test obtains the exciting force spectrum on the end effector mounting bracket in the working process of palletizing robot. Main experimental equipment still includes LMS TEST.LAB system, LMS SC305 system, three one-way piezoelectric acceleration sensors, and computers. Because vibration is passed on to the end effector mounting bracket mainly through bolts, therefore the test points, a total of 3, were set near the bolts on the end effector. Measurement: the test points are set on the center of three bearing caps, as shown in Figure6. Each test point is measured at the acceleration signal of x, y, and z directions. Data acquisition parameters setting: sampling time is 30 s, frequency bandwidth is 1024 Hz, and resolution is 0.125 Hz. The acceleration spectrum in the +x direction of #1 test point, for example, is shown in Figure 7.Figure 6
Test points location.Figure 7
Acceleration spectrum in the +x direction of #1 test point.


### 2.5. Frequency Response Analysis of the End Effector Mounting Bracket

In accordance with the position of test points, the measured acceleration spectrums are applied on the corresponding position points of the finite element model, respectively, along thex, y, and z directions, and then the frequency response results are solved. By comparing and analyzing maximum vibration velocity (4.08 mm/s) at node A in the red zone on the model, as shown in Figure 8, this may cause fatigue damage of parts. The frequency response curve of node A, as shown in Figure 9, shows that there are two peak vibration velocities: #1 and #2; the #1 peak vibration velocity is the maximum (4.08 mm/s) that occurs at a frequency of about 430 Hz, in the vicinity of the first natural frequency of the end effector mounting bracket. In order to reduce the impact of vibration on the reliability of the parts and motion accuracy of the robot, choose the maximizing first natural frequency as one of the optimization goals.Figure 8
Frequency response of initial model.Figure 9
Frequency response curve of nodeA.


### 2.6. Static Analysis of the End Effector Mounting Bracket

In palletizing process, the end effector carrying load undergoes acceleration and deceleration process; therefore, inertia force should be taken into account in the static analysis. In this paper, palletizing robot picks material from the conveyor belt and places it on the specified position on the tray along a plane “door” shaped motion trajectory. In order to avoid shock and vibration phenomenon caused by the sudden change of velocity and shorten the operation cycle to increase efficiency, two turning points (Pm and Pn) are replaced by two arcs, with motion trajectory of the end effector as shown in Figure 10. The end effector picks the material from P1, followed by P2, P3, P4, P5, and P6, and places it on the specified location tray. The establishment of coordinate system o-xyz at the P1 point (y direction perpendicular to the paper surface outward) and the coordinate value of key points on the target trajectory is as follows (unit: mm): P1=[0,0,0], Pm=[0,0,700], Pn=[-2400,0,700], and P6=[-2400,0,0].Figure 10
Motion trajectory of the palletizing robot end diagram.The modified trapezoid motion rule is adopted [1]; its distance, velocity, and acceleration have varying tendency with time as shown in Figure 11. During operation, maximum acceleration amax is 3.2 m/s2 and maximum velocity vmax is 1.2 m/s.Figure 11
Modified trapezoidal rule.Maximum acceleration, which occurs in vertically accelerated lifting and decelerated falling process of the materials, will make the material in the overweight state. This moment, the end effector mounting bracket bears maximum force in the process of the whole movement; therefore, this maximum force is used as load boundary condition of static analysis, and its value can be calculated as(1)Fmax=Mg+amax,where Fmax is the maximum force (N);  M is the total mass including maximum load mass and the total mass of other parts of the end effector without end effector mounting bracket, with a value of 132.74 kg; g is the gravitational acceleration, with a value of 9.8 m/s2.By (1), the maximum load on the motor installation surface of the end effector mounting bracket size is 1725.62 N, and its direction is the gravitational direction. The displacement constraints are cylindrical surface type constraints in three bearing holes.The stress nephogram and deformation nephogram are obtained by the statics analysis, as shown in Figures12 and 13, respectively. The maximum stress is about 6.8 MPa, far less than the allowable stress values of material; the maximum deformation is 0.042 mm, because the greater the deformation, the lower the stiffness of structure, and therefore the stiffness of structure can be characterized by the size of deformation [30]. Statics analysis shows that the end effector mounting bracket is of great potential in structure optimization.Figure 12
Stress nephogram.Figure 13
Deformation nephogram.


### 2.7. Establishment of Multicriteria Optimization Model for the End Effector Mounting Bracket

#### 2.7.1. The Choice of Design Variables

(1) Initial Choice of Design Variables. Initially, 9 structure parameters X=(x1,x2,…,x9) of the end effector mounting bracket were selected as shown in Figure 14; their names, initial values, and value range are shown in Table 5. Structural parameters are all local parameters, independent of each other, without changing the overall size and the assembly dimensions with other parts.Table 5
The initial value and range of design variables, mm.
Design variable Name Initial value Value range x 1 Small hole diameter of high upright plate 50 50~75 x 2 Large hole diameter of high upright plate 80 80~120 x 3 Arc upright plate thickness 15 6~15 x 4 Bottom plate thickness 25 8~25 x 5 Posterior slope plate thickness 15 6~15 x 6 Elliptical hole major axis of middle plate 100 100~180 x 7 Elliptical hole minor axis of middle plate 50 50~80 x 8 Middle plate thickness 15 6~15 x 9 Stiffener thickness of high upright plate 15 6~15Figure 14
Initial choice of design variables.(2) Sensitivity Analysis. In the process of optimizing design, if the design variables are too many, the key structure parameters can be neglected, resulting in the pseudo optimization results, and the calculated quantity is large and the efficiency is low. Through the sensitivity analysis, some structural design parameters, which have significantly influenced the structure response parameter, are selected as the design variable from the initial selection [30].There are many kinds of sensitivity analysis methods; among them, sensitivity analysis method based on experimental data combined with Spearman rank correlation coefficient method is a reliable and effective method [31, 32].(i) Design of Experiments (DOE) Method. Experimental data can be obtained by DOE method such as Central Composite Design (CCD), Doehlert Matrix (DM), and the Box-Behnken design (BBD) method [33]. Compared with the CCD and DM methods, Box-Behnken design method is more efficient and economical; therefore, Box-Behnken method [34] is used to estimate coefficient of the second-order polynomial response surface model. The Box-Behnken design method is an incomplete three-level part factorial experiment design method, is of nearly rotatable symmetry, is composed of multiple orthogonal cubes, and includes a central point. Schematic diagram of the three-level Box-Behnken experimental design plan is shown in Figure 15. The advantage of the Box-Behnken design is that it can avoid emergence extreme point [35].Figure 15
Schematic diagram of Box-Behnken experimental design plan.(ii) Sensitivity Analysis Based on DOE Combined with Spearman Rank Correlation Coefficient Method. Suppose sensitivity of structural response y (such as maximum stress) to structural design parameters x1,x2,x3,…,xm is examined. According to the initial value and value range of structural design parameters, we use the design of experiments (DOE) method to obtain a set of sample points of these structural design parameters x1i,x2i,x3i,…,xmi (i=1,2,…,n), and n response values of maximum stress y1,y2,y3,…,yn are obtained by n times calculating simulations. m data pairs consist of the sample values xj1,xj2,xj3,…,xjn that are calculated n times by the jth (j=1,2,3,…,m) structural design parameters, and the response values of maximum stress y1,y2,y3,…,yn by calculating simulations n times are as follows:(2)xj1y1,xj2y2,…,xjnynj=1,2,…,m.The correlation between thejth structural design parameters and the maximum stress is investigated by the Spearman rank correlation coefficient method. Respectively, the n sample values xj1,xj2,xj3,…,xjn of the jth structural design parameters and the response value y1,y2,y3,…,yn of the maximum stress are sorted in descending order, respectively, using Rji,Qi to indicate the position number of the original data xji, yi (i=1,2,…,n) in the sorted list, where Rji,Qi are called the rank of xji, yi A set of new data pairs are formed by Rji, Qi instead of xji, yi as follows: (3)Rj1Q1,Rj2Q2,…,RjnQn.Then Spearman rank correlation coefficient of the original data pairs is(4)rsj=∑i=1nRji-R-Qi-Q-∑i=1nRji-R-2∑i=1nQi-Q-2,j=1,2,3,…,m,where (5)R-=∑i=1nRjin=n+12,Q-=∑i=1nQin=n+12.Simplification of formula (1):(6)rsj=1-6∑i=1nRji-Qi2nn2-1.Rank correlation coefficientrsj, in the range of [-1,1], reflects the correlation between the maximum stress response y and the jth structural parameters. That rsj has a positive value indicates a monotonically increasing relationship between the two; that rsj=1 indicates perfect positive correlation; that rsj has a negative value indicates a monotonically decreasing relationship between the two; that rsj=-1 indicates perfect negative correlation; the greater the absolute value rsj, the greater the correlation, that is, the more sensitive the structure response y to the jth structural parameters; rsj=0 indicates that there is no relationship between two variables.By that analogy, the sensitivity of structural response parameters to all selected structural design parameters can be obtained.According to the value range of initial design variables as well as Box-Behnken design method generating experiment design matrix table, with a total of 97 sets, the mass (m, unit: kg), the first natural frequency (f1, unit: Hz), the maximum stress (σmax, unit: MPa), and the maximum deformation (δmax, unit: mm) are calculated, respectively. The experiment design matrix and corresponding results are obtained as shown in Tables 6 and 7.Table 6
Experiment design matrix, mm.
Number Design variables x 1 x 2 x 3 x 4 x 5 x 6 x 7 x 8 x 9 1 62.5 100 10.5 16.5 10.5 140 65 10.5 10.5 2 50 100 10.5 8 10.5 140 50 10.5 10.5 3 75 100 10.5 8 10.5 140 50 10.5 10.5 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ 97 62.5 100 15 16.5 15 140 80 10.5 10.5Table 7
Experiment design results.
Number Output result m f 1 δ m a x σ m a x 1 7.5701 391.18 0.0603 7.7425 2 7.0498 399.94 0.0765 9.6833 3 6.9479 399.60 0.0773 9.5524 ⋯ ⋯ ⋯ ⋯ ⋯ 97 8.0109 396.10 0.0499 7.4748Then, the Spearman rank correlation coefficient method is used to obtain the mass sensitivity, the first natural frequency sensitivity, the maximum stress sensitivity, and the maximum deformation sensitivity of the end effector mounting bracket, as shown in Figure16.Figure 16
Sensitivity analysis of structural design parameters to response parameters.
(a)
Mass sensitivity (b)
The maximum deformation sensitivity (c)
The maximum stress sensitivity (d)
The first natural frequency sensitivityHere, five structural design parameters, whose every kind of sensitivity value is greater than 0.2, are selected as design variables. The design variables are identified as(7)X=x2,x3,x4,x6,x8.


#### 2.7.2. Constraint Condition

Constraint conditions are values range of every design variable.(1) Design Variable x 2 and x 6 Are the Diameter of the Hole. The greater the hole, the lighter the mass, therefore initial values of these holes diameter are set to the minimum value, their maximum values are less than boundary dimensions of structure.(2) Design Variables x 3 , x 4,andx8Are the Thickness of Plates. The smaller the thickness, the lighter the mass; therefore, the initial values of these plates’ thicknesses are set to the maximum value, and their minimum value is set to the minimum thickness value of cast aluminum alloy part [26], as shown in Table 5.


#### 2.7.3. The Establishment of Objective Function

(1) Response Surface Methodology (RSM). The end effector mounting bracket is an irregular part; by the structural parameters, directly constructing accurate objective functions is difficult; therefore, Response Surface Methodology (RSM) [36] is introduced, using numerical approximation and functional fit method to construct objective functions. The basic principle of RSM is to use DOE [37] method experiments with design variables to obtain response surface model of objective functions and predict the response value of nonexperimental points. Suppose a function of design variables X and response variables F(X) can be described as(8)FX=F~X+ε=∑i=0LβiϕiX+ε,where F(X) is the actual response function of the optimization objective for the design variables, F~(X) is the approximate response function, namely, response surface model, L is the number of basis functions, ϕi(X) is a polynomial function of design variable X, and ε is the random error of the approximate value and the actual value.The approximate models of objective functions are constructed by the second-order polynomial response surface model as shown in(9)F~X=β0+∑i=1nβixi+∑i=1nβiixi2+∑i<jnβijxixj,where xi are input variables, β0 is the constant, and βi is a variable coefficient of regression model.Using the least square method solves the coefficient of (9): (10)A=β0,β1,…,βLT.Then,(11)A=XTX-1XTY,where X is the basis function matrix, which is described as (12)X=ϕ1X1ϕ2X1⋯ϕLX1ϕ1X2ϕ2X2⋯ϕLX2⋯⋯⋱⋯ϕ1XQϕ2XQ⋯ϕLXQ.(2) The Establishment of Objective Functions Based on the RSM. With minimal mass, minimizing the maximal stress, minimizing the maximal deformation, and maximizing the first natural frequency as optimization objectives, therefore objective functions are established as (13)FmX=minX∈DmX,Ff1X=maxX∈Df1X,FσX=minX∈DσmaxX,FδX=minX∈DδmaxX.According to the value range of selected design variables as well as the Box-Behnken design method generating experiment design matrix table, with a total of 41 sets, the mass (m, unit: kg); the first natural frequency (f1, unit: Hz), the maximum stress (σmax, unit: MPa), and maximum deformation (δmax, unit: mm) are calculated, respectively. The experiment design matrix and its results are obtained as shown in Tables 8 and 9, respectively.Table 8
Experiment design matrix.
Number Design variable (mm) x 2 x 3 x 4 x 6 x 8 1 120 15 16.5 140 10.5 2 120 6 16.5 140 10.5 3 80 15 16.5 140 10.5 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ 41 100 10.5 16.5 140 10.5Table 9
Experiment design results.
Number Output results m f 1 δ m a x σ m a x 1 8.0534 402.72 0.0508 5.9895 2 7.3791 399.45 0.0774 7.1030 3 8.3144 398.43 0.0489 5.9539 ⋯ ⋯ ⋯ ⋯ ⋯ 41 7.8617 399.85 0.0587 6.0313With the design variables in Table8 as inputs and the calculated results in Table 9 as outputs, the response surface models of the four objective functions are obtained as in (14)~(17), respectively, where m(X) is the mass function, f1(X) is the first natural frequency function, σmax(X) is the maximum stress function, and δmax(X) is the maximum deformation function. (14)mX=5.56-1.2×10-5x2+0.1x3+0.08x4-6.88×10-5x6+0.06x8-3.26×10-5x22-9.36×10-5x32+3.06×10-7x42+2.44×10-7x62+1.04×10-5x82+1.15×10-18x2x3+1.47×10-8x2x4+3.13×10-9x2x6+2.78×10-8x2x8-0.001x3x4+1.39×10-8x3x6+1.23×10-7x3x8+2.18×10-19x4x6+9×10-18x4x8-0.0001x6x8,(15)f1X=345.91+0.4x2+1.84x3+0.18x4+0.05x6+3.49x8-0.0004x22-0.09x32-0.01x42-0.001x62-0.03x82-0.008x2x3-0.005x2x4-0.0002x2x6+0.002x2x8+0.03x3x4+0.01x3x6-0.07x3x8-0.0013x4x6+0.007x4x8-0.002x6x8,(16)δmaxX=0.18+2.011x2-0.01x3-0.006x4+6.18×10-5x6-0.0007x8+3.08×10-7x22+0.0002x32+6.76×10-5x42+8.37×10-8x62-2.98×10-7x82-1.92×10-6x2x3-2.5×10-7x2x4+9.37×10-9x2x6-3.89×10-7x2x8+0.0002x3x4-5.58×10-6x3x6+4.4×10-5x3x8+9.26×10-7x4x6-6.34×10-6x4x8-1.39×10-8x6x8,(17)σmaxX=13.1+0.008x2-0.85x3-0.33x4+0.004x6+0.09x8+4.58×10-5x22+0.023x32+0.004x42+5.6×10-5x62+0.005x82-0.0001x2x3-0.0001x2x4-2.19x2x6-0.0009x2x8+0.016x3x4-0.0004x3x6+0.0028x3x8+8.39×10-5x4x6-0.002x4x8-0.001x6x8.In order to validate the credibility of the response surface models, reselect 10 sets’ value of design parameters in their value range by the Box-Behnken method and calculate; the calculated results were compared with response surface models predicted results as shown in Table10. Due to the fact that the mass error is 0, therefore it is not listed in the table. Table 10 shows all the errors are within 2%; the credibility of the above response surface models meets design requirements, which can be used as objective functions for the optimization design [38].Table 10
Comparison between the calculation results and predicted results of the approximate model.
Number Calculated results Response surface model predicted results Relative error (%) f 1 δ m a x σ m a x f 1 δ m a x σ m a x f 1 δ m a x σ m a x 1 398.73 0.0504 6.0405 399.05 0.0497 5.9992 0.0816 1.4429 0.6840 2 382.31 0.0552 6.2636 382.54 0.0553 6.3163 0.0591 0.1819 0.8411 3 388.72 0.0534 6.0975 389.01 0.0532 5.9920 0.0722 0.4544 1.7298 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ 10 395.56 0.0576 6.0904 395.89 0.0570 6.0237 0.0833 1.0272 1.0945


#### 2.7.4. Determine the Weight Coefficients

This paper presents an analytic hierarchy process based on finite element analysis results (FEA + AHP) method to determine the weight coefficients of each optimization objective.(1) Analytic Hierarchy Process (AHP). The basic principle and steps of the AHP are as follows [39].(i) Establish Hierarchical Structure Model. Generally, a hierarchical structure model includes three levels: goal level, criteria level, and policy level, as shown in Figure 17. Hierarchy number is unrestricted; generally, the number of the elements is not more than nine in each level.Figure 17
Hierarchical structure model.(ii) Construct Pairwise Comparison Matrix. If a level has n factors, then compare their importance to a certain factor of the above level. The 1–9 scale method is used to obtain mark results as shown in Table 11. The pairwise comparison matrix A is obtained according to the mark results, with matrix A as shown in (18)A=aijn×n=a11a12⋯a1na21a22⋯a2n⋯⋯⋯⋯an1an2⋯ann.Table 11
Criteria of 1~9 scale.
Scale Scale ofi compared to j 1 Equally important 3 Weakly important 5 Strongly important 7 Very strongly important 9 Extremely important 2, 4, 6, 8 Intermediate values betweentwo adjacency values Reciprocal If marki to j is aij, then mark j to i is aij=1/aij(iii) Hierarchy Sorting and Consistency Checking. (a) Calculate the consistency index (CI): (19)CI=λmax-nn-1,where λmax is the maximum eigenvalue of the comparison matrix.(b) The random consistency index RI is shown in Table12.Table 12
Random index.
n 1 2 3 4 ⋯ 10 RI 0 0 0.58 0.90 ⋯ 1.51(c) Calculate the consistency ratio CR:(20)CR=CIRI.WhenCR<0.1, the consistency of the comparison matrix is acceptable; otherwise, the comparison matrix should be properly modified until CR<0.1.(iv) Solving the Maximum Eigenvalue λ m a x Corresponding Eigenvector W ∗ = ( ω 1 , ω 2 , … , ω n ) T of the Comparison Matrix Which Passes the Consistency Checking. And then the weight vector W is obtained by the normalization to the eigenvector W∗.(2) Determine the Weight Coefficients(i) Establishing the Hierarchical Structure Model Shown in Figure18. The goal level is structure optimization, using A to represent it. The criteria level is mainly composed of the mass, the maximum deformation, the first natural frequency, and the maximum stress, respectively, using B1, B2, B3, and B4 to represent them. The policy level is mainly composed of five design variables, respectively, using C1, C2, C3, C4, and C5 to represent them.Figure 18
Hierarchical structure model.(ii) Construct the Comparison Matrix. Based on the results of finite element analysis, the importance of each performance index is investigated. The main objective of this study is the lightweight of the end effector mounting bracket; therefore, reducing mass is the main goal. The importance of the first natural frequency, deformation, and stress is discussed. It is known from the frequency response analysis that vibration frequency corresponds to the maximum vibration velocity located in the vicinity of the first natural frequency, easily causing vibration fatigue of parts; consequently, its importance is ranked the second; as can be seen from the statics analysis, the maximum deformation is very small; however, in order to reduce the impact of the deformation on the positioning accuracy of the palletizer robot, the maximum deformation should be smaller to be better; therefore, the importance of deformation is ranked the third; the maximum stress is far less than the allowable stress; therefore, its importance comes in the last row. Accordingly, the pairwise comparison results of each factor in the criterion layer are shown in Table 13. And get the comparison matrix as shown in(21)A=1537151122132131712131.Table 13
Pairwise comparison results of each factor.
A B 1 B 2 B 3 B 4 B 1 1 5 3 7 B 2 1/5 1 1/2 2 B 3 1/3 2 1 3 B 4 1/7 1/2 1/3 1The eigenvalue of maximum of matrixA is (22)λmax=4.0192;λ m a x corresponding eigenvector is (23)W∗=0.9142,0.1912,0.3392,0.1124T.(iii) Consistency Checking Consistency index:CI=0.0064<0.1. Random consistency index:RI=0.9. Consistency ratio:CR=0.007<0.1.Accordingly, it can be determined whether comparison matrix A passes the consistency checking. Then, we carry out the eigenvector normalization process to obtain the weight vector of minimum mass, the maximum deformation, the maximum first natural frequency, and maximum stress:(24)W=0.5872,0.1228,0.2179,0.0722T.


### 2.8. Solving Multicriteria Optimization Problem

#### 2.8.1. Elitist Nondominated Sorting Genetic Algorithm (NSGA-II)

NSGA-II algorithm is an improved algorithm of NSGA, which was proposed by Deb et al. [40]. NSGA-II algorithm uses fast nondominated sorting procedure, elite strategy, and nonparameter ecological niche operator, overcoming the shortcomings of the traditional NSGA such as high computational complexity, nonelite strategy, and the need for specially designated shared radius. NSGA-II algorithm has been widely used due to its good performance in exploring; Pareto set, which is obtained by NSGA-II algorithm, has good accuracy and dispersity. The NSGA-II procedure is shown in Figure 19 [7, 40].Figure 19
The main flow chart of NSGA-II.Step 1.
Initially, a random parent populationP0 is created. The population is sorted based on the nondomination. Each solution is assigned a fitness (or rank) to its nondomination level; binary tournament selection, recombination, and mutation operators are used to create a child population Q0 of size N and set t=0.Step 2.
A combined populationRt=Pt∪Qt is formed. Here, Pt is the parent population in tth iteration; Qt is the child population in tth iteration; t denotes the iteration or generation number; the population Rt is of size 2N. Then, Rt is sorted according to nondomination. Since all previous and current population members are included in Rt, elitism is ensured. And then nondominated fronts F1,F2,…,Fi are obtained.Step 3.
SortFi using the crowded-comparison operator <n in descending order and choose the best solutions to form a new population Pt+1.Step 4.
The new populationPt+1 of size N is used for selection, crossover, and mutation to create a new population Qt+1 of size N.Step 5.
If the termination conditions are satisfied, it is over; otherwise,t=t+1; go to Step 2.Parameter configuration of NSGA-II algorithm: population size is 40; number of generations is 200; crossover rate is 0.9; cross distribution index is 10; mutation distribution index is 20.Justification of NSGA-II algorithm parameters selection: firstly, population size directly affects astringency or computational efficiency. If population size is too small, optimization results are likely to converge to local optimal solution. If population size is too large, the calculations will take a long time. Therefore, population size usually is selected in the range 20~200. Secondly, the meaning of number of generations is that NSGA-II algorithm will be terminated when the number of evolutionary generations reaches this set value. Generally, the number of generations is selected in the range 100~500 [41]. In this work, the above two parameters are determined by some tests for permutation and combination of different values of population size and number of generations in their respective ranges. The comparison test results indicate that there is only small difference in the decimal part among the group results; however, precision of results is improved with the increase of population size and number of generations, and the calculation time has become longer. Population size = 40 and number of generations = 200 are determined by taking into account the calculation precision and time. Thirdly, crossover rate controls frequentness of crossover operation. If crossover rate is too high, generation gap is likely to happen. However, if crossover rate is too low, many individuals will be prone to duplicating in the next generation directly, leading to stagnation during searching. Therefore, crossover rate is usually selected in the range 0.4~0.99 [41]. In some literatures, crossover rate is 0.9 [38, 42]. Finally, cross distribution index and mutation distribution index both are nonnegative number and user-defined rule, and the higher their values, the smaller the probability of offspring individuals far away from parent individuals, and vice versa. Therefore, the selected cross distribution index and mutation distribution index, respectively, are 10 and 20; they are smaller, in order to enhance searching ability of algorithm [38, 43].


## 2.1. Establishment of Finite Element Model of the End Effector Mounting Bracket

Before the finite element model is established, the 3D model is simplified by removing chamfer, fillet, screw holes, and others of the part. The simplified model is shown in Figure2. The material property and mesh setting are shown in Table 1 [26] and Table 2, respectively. The finite element model is shown in Figure 3.Table 1
Material property.
Item Property Material ZL110 Elastic modulus (N/m2) 6.89 × 1010 Poisson ratio 0.35 Shear modulus (N/m2) 3.19 × 108 Density (kg/m3) 2770 Tensile strength (N/m2) 1.65 × 108Table 2
Mesh setting.
Item Property Mesh type Tetrahedral Minimum element size (mm) 4.96532 Maximum element size (mm) 14.8961 Total number of nodes 19995 Total number of units 10884Figure 2
Simplified model of end effector mounting bracket.Figure 3
Finite element model of end effector mounting bracket.The following is a justification of the above mesh setting. First of all, element selection: tetrahedral mesh has been gaining popularity for solving finite element analysis with its simple, flexible, and strong adaptability to complex boundary (Chunge et al. [27]). The end effector mounting bracket is an irregular and complicated part; therefore, tetrahedral mesh is selected in this work. Second, number of elements: in order to ensure the accuracy of finite element analysis, the highest quality elements are adopted in the finite element software; therefore, the number of elements, total number of nodes, minimum element size, and maximum element size are obtained automatically.


## 2.2. Finite Element Model Accuracy Verification Based on the Modal Test

The modal test method is used to verify the accuracy of the finite element model. Experimental equipment includes LMS TEST.LAB vibration and noise analysis system, LMS SC305 type vibration and noise data acquisition system, three Lance LC0152T type unidirectional piezoelectric acceleration sensors, B&K Impact Hammer Type 8207, signal cable, IBM computer, crane, and elastic rope. We use single-point excitation and multipoint response test method that selects 21 pick points and 1 knock point (O point) on the part and then use impact hammer knock on the O point, respectively, along -x, +y, and +z directions; meanwhile, x, y, and z direction responses of every pick point (a total of 63 sets of data) are measured with three acceleration sensors, and measurement signals are acquired by the LMS SC305 system and are processed by the LMS TEST.LAB system; the test modal of the end effector mounting bracket is obtained afterwards. To measure the free modal of the part, it is overhung by a greatly flexible elastic rope to make it close to the free boundary conditions. The experimental environment is shown in Figure 4.Figure 4
The experimental environment.The free modal of the end effector mounting bracket is calculated by the finite element method. Comparing the FEM model with test modal for the first three modal shapes and the first five modal frequencies is shown in Figure5 and Table 3, respectively, where Ftestmodal is test modal frequency and Ffemmodal is finite element modal frequency.Table 3
Comparison between the test model and the FEM model for the first five modal frequencies (Hz).
Modal order F t e s t   m o d a l F f e m   m o d a l Relative error% 1 229.11 228.34 0.3404 2 448.41 449.75 0.3011 3 757.03 732.32 3.2628 4 789.54 768.58 2.6547 5 932.24 942.74 1.1263Figure 5
Comparison between the test model and FEM model for the first three modal shapes.
(a)
Comparison between the test model and FEM model for the first modal shape (b)
Comparison between the test model and FEM model for the second modal shape (c)
Comparison between the test model and FEM model for the third modal shapeThe comparison results indicate that the first three modal shapes of the calculating modal and test modal are consistent, and the relative error between the first five modal frequencies of the calculating modal and test modal is less than 5%. This shows that the accuracy of the established finite element model meets the requirements, and it can be used for the subsequent calculating simulation [28, 29].


## 2.3. Constraint Modal Analysis of the End Effector Mounting Bracket

According to the assembly relation of end effector mounting bracket and other parts, cylindrical surface constraints are applied in three bearing holes of its finite element model, and the first four modals are calculated by making use of Block Lanczos method. The first four natural frequencies are shown in Table4.Table 4
The first four natural frequencies of the end effector mounting bracket, Hz.
Order number 1 2 3 4 Frequency 409.57 802.07 1021.2 1312.1


## 2.4. Vibration Response Test

Using the vibration response test obtains the exciting force spectrum on the end effector mounting bracket in the working process of palletizing robot. Main experimental equipment still includes LMS TEST.LAB system, LMS SC305 system, three one-way piezoelectric acceleration sensors, and computers. Because vibration is passed on to the end effector mounting bracket mainly through bolts, therefore the test points, a total of 3, were set near the bolts on the end effector. Measurement: the test points are set on the center of three bearing caps, as shown in Figure6. Each test point is measured at the acceleration signal of x, y, and z directions. Data acquisition parameters setting: sampling time is 30 s, frequency bandwidth is 1024 Hz, and resolution is 0.125 Hz. The acceleration spectrum in the +x direction of #1 test point, for example, is shown in Figure 7.Figure 6
Test points location.Figure 7
Acceleration spectrum in the +x direction of #1 test point.


## 2.5. Frequency Response Analysis of the End Effector Mounting Bracket

In accordance with the position of test points, the measured acceleration spectrums are applied on the corresponding position points of the finite element model, respectively, along thex, y, and z directions, and then the frequency response results are solved. By comparing and analyzing maximum vibration velocity (4.08 mm/s) at node A in the red zone on the model, as shown in Figure 8, this may cause fatigue damage of parts. The frequency response curve of node A, as shown in Figure 9, shows that there are two peak vibration velocities: #1 and #2; the #1 peak vibration velocity is the maximum (4.08 mm/s) that occurs at a frequency of about 430 Hz, in the vicinity of the first natural frequency of the end effector mounting bracket. In order to reduce the impact of vibration on the reliability of the parts and motion accuracy of the robot, choose the maximizing first natural frequency as one of the optimization goals.Figure 8
Frequency response of initial model.Figure 9
Frequency response curve of nodeA.


## 2.6. Static Analysis of the End Effector Mounting Bracket

In palletizing process, the end effector carrying load undergoes acceleration and deceleration process; therefore, inertia force should be taken into account in the static analysis. In this paper, palletizing robot picks material from the conveyor belt and places it on the specified position on the tray along a plane “door” shaped motion trajectory. In order to avoid shock and vibration phenomenon caused by the sudden change of velocity and shorten the operation cycle to increase efficiency, two turning points (Pm and Pn) are replaced by two arcs, with motion trajectory of the end effector as shown in Figure 10. The end effector picks the material from P1, followed by P2, P3, P4, P5, and P6, and places it on the specified location tray. The establishment of coordinate system o-xyz at the P1 point (y direction perpendicular to the paper surface outward) and the coordinate value of key points on the target trajectory is as follows (unit: mm): P1=[0,0,0], Pm=[0,0,700], Pn=[-2400,0,700], and P6=[-2400,0,0].Figure 10
Motion trajectory of the palletizing robot end diagram.The modified trapezoid motion rule is adopted [1]; its distance, velocity, and acceleration have varying tendency with time as shown in Figure 11. During operation, maximum acceleration amax is 3.2 m/s2 and maximum velocity vmax is 1.2 m/s.Figure 11
Modified trapezoidal rule.Maximum acceleration, which occurs in vertically accelerated lifting and decelerated falling process of the materials, will make the material in the overweight state. This moment, the end effector mounting bracket bears maximum force in the process of the whole movement; therefore, this maximum force is used as load boundary condition of static analysis, and its value can be calculated as(1)Fmax=Mg+amax,where Fmax is the maximum force (N);  M is the total mass including maximum load mass and the total mass of other parts of the end effector without end effector mounting bracket, with a value of 132.74 kg; g is the gravitational acceleration, with a value of 9.8 m/s2.By (1), the maximum load on the motor installation surface of the end effector mounting bracket size is 1725.62 N, and its direction is the gravitational direction. The displacement constraints are cylindrical surface type constraints in three bearing holes.The stress nephogram and deformation nephogram are obtained by the statics analysis, as shown in Figures12 and 13, respectively. The maximum stress is about 6.8 MPa, far less than the allowable stress values of material; the maximum deformation is 0.042 mm, because the greater the deformation, the lower the stiffness of structure, and therefore the stiffness of structure can be characterized by the size of deformation [30]. Statics analysis shows that the end effector mounting bracket is of great potential in structure optimization.Figure 12
Stress nephogram.Figure 13
Deformation nephogram.


## 2.7. Establishment of Multicriteria Optimization Model for the End Effector Mounting Bracket

### 2.7.1. The Choice of Design Variables

(1) Initial Choice of Design Variables. Initially, 9 structure parameters X=(x1,x2,…,x9) of the end effector mounting bracket were selected as shown in Figure 14; their names, initial values, and value range are shown in Table 5. Structural parameters are all local parameters, independent of each other, without changing the overall size and the assembly dimensions with other parts.Table 5
The initial value and range of design variables, mm.
Design variable Name Initial value Value range x 1 Small hole diameter of high upright plate 50 50~75 x 2 Large hole diameter of high upright plate 80 80~120 x 3 Arc upright plate thickness 15 6~15 x 4 Bottom plate thickness 25 8~25 x 5 Posterior slope plate thickness 15 6~15 x 6 Elliptical hole major axis of middle plate 100 100~180 x 7 Elliptical hole minor axis of middle plate 50 50~80 x 8 Middle plate thickness 15 6~15 x 9 Stiffener thickness of high upright plate 15 6~15Figure 14
Initial choice of design variables.(2) Sensitivity Analysis. In the process of optimizing design, if the design variables are too many, the key structure parameters can be neglected, resulting in the pseudo optimization results, and the calculated quantity is large and the efficiency is low. Through the sensitivity analysis, some structural design parameters, which have significantly influenced the structure response parameter, are selected as the design variable from the initial selection [30].There are many kinds of sensitivity analysis methods; among them, sensitivity analysis method based on experimental data combined with Spearman rank correlation coefficient method is a reliable and effective method [31, 32].(i) Design of Experiments (DOE) Method. Experimental data can be obtained by DOE method such as Central Composite Design (CCD), Doehlert Matrix (DM), and the Box-Behnken design (BBD) method [33]. Compared with the CCD and DM methods, Box-Behnken design method is more efficient and economical; therefore, Box-Behnken method [34] is used to estimate coefficient of the second-order polynomial response surface model. The Box-Behnken design method is an incomplete three-level part factorial experiment design method, is of nearly rotatable symmetry, is composed of multiple orthogonal cubes, and includes a central point. Schematic diagram of the three-level Box-Behnken experimental design plan is shown in Figure 15. The advantage of the Box-Behnken design is that it can avoid emergence extreme point [35].Figure 15
Schematic diagram of Box-Behnken experimental design plan.(ii) Sensitivity Analysis Based on DOE Combined with Spearman Rank Correlation Coefficient Method. Suppose sensitivity of structural response y (such as maximum stress) to structural design parameters x1,x2,x3,…,xm is examined. According to the initial value and value range of structural design parameters, we use the design of experiments (DOE) method to obtain a set of sample points of these structural design parameters x1i,x2i,x3i,…,xmi (i=1,2,…,n), and n response values of maximum stress y1,y2,y3,…,yn are obtained by n times calculating simulations. m data pairs consist of the sample values xj1,xj2,xj3,…,xjn that are calculated n times by the jth (j=1,2,3,…,m) structural design parameters, and the response values of maximum stress y1,y2,y3,…,yn by calculating simulations n times are as follows:(2)xj1y1,xj2y2,…,xjnynj=1,2,…,m.The correlation between thejth structural design parameters and the maximum stress is investigated by the Spearman rank correlation coefficient method. Respectively, the n sample values xj1,xj2,xj3,…,xjn of the jth structural design parameters and the response value y1,y2,y3,…,yn of the maximum stress are sorted in descending order, respectively, using Rji,Qi to indicate the position number of the original data xji, yi (i=1,2,…,n) in the sorted list, where Rji,Qi are called the rank of xji, yi A set of new data pairs are formed by Rji, Qi instead of xji, yi as follows: (3)Rj1Q1,Rj2Q2,…,RjnQn.Then Spearman rank correlation coefficient of the original data pairs is(4)rsj=∑i=1nRji-R-Qi-Q-∑i=1nRji-R-2∑i=1nQi-Q-2,j=1,2,3,…,m,where (5)R-=∑i=1nRjin=n+12,Q-=∑i=1nQin=n+12.Simplification of formula (1):(6)rsj=1-6∑i=1nRji-Qi2nn2-1.Rank correlation coefficientrsj, in the range of [-1,1], reflects the correlation between the maximum stress response y and the jth structural parameters. That rsj has a positive value indicates a monotonically increasing relationship between the two; that rsj=1 indicates perfect positive correlation; that rsj has a negative value indicates a monotonically decreasing relationship between the two; that rsj=-1 indicates perfect negative correlation; the greater the absolute value rsj, the greater the correlation, that is, the more sensitive the structure response y to the jth structural parameters; rsj=0 indicates that there is no relationship between two variables.By that analogy, the sensitivity of structural response parameters to all selected structural design parameters can be obtained.According to the value range of initial design variables as well as Box-Behnken design method generating experiment design matrix table, with a total of 97 sets, the mass (m, unit: kg), the first natural frequency (f1, unit: Hz), the maximum stress (σmax, unit: MPa), and the maximum deformation (δmax, unit: mm) are calculated, respectively. The experiment design matrix and corresponding results are obtained as shown in Tables 6 and 7.Table 6
Experiment design matrix, mm.
Number Design variables x 1 x 2 x 3 x 4 x 5 x 6 x 7 x 8 x 9 1 62.5 100 10.5 16.5 10.5 140 65 10.5 10.5 2 50 100 10.5 8 10.5 140 50 10.5 10.5 3 75 100 10.5 8 10.5 140 50 10.5 10.5 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ 97 62.5 100 15 16.5 15 140 80 10.5 10.5Table 7
Experiment design results.
Number Output result m f 1 δ m a x σ m a x 1 7.5701 391.18 0.0603 7.7425 2 7.0498 399.94 0.0765 9.6833 3 6.9479 399.60 0.0773 9.5524 ⋯ ⋯ ⋯ ⋯ ⋯ 97 8.0109 396.10 0.0499 7.4748Then, the Spearman rank correlation coefficient method is used to obtain the mass sensitivity, the first natural frequency sensitivity, the maximum stress sensitivity, and the maximum deformation sensitivity of the end effector mounting bracket, as shown in Figure16.Figure 16
Sensitivity analysis of structural design parameters to response parameters.
(a)
Mass sensitivity (b)
The maximum deformation sensitivity (c)
The maximum stress sensitivity (d)
The first natural frequency sensitivityHere, five structural design parameters, whose every kind of sensitivity value is greater than 0.2, are selected as design variables. The design variables are identified as(7)X=x2,x3,x4,x6,x8.


### 2.7.2. Constraint Condition

Constraint conditions are values range of every design variable.(1) Design Variable x 2 and x 6 Are the Diameter of the Hole. The greater the hole, the lighter the mass, therefore initial values of these holes diameter are set to the minimum value, their maximum values are less than boundary dimensions of structure.(2) Design Variables x 3 , x 4,andx8Are the Thickness of Plates. The smaller the thickness, the lighter the mass; therefore, the initial values of these plates’ thicknesses are set to the maximum value, and their minimum value is set to the minimum thickness value of cast aluminum alloy part [26], as shown in Table 5.


### 2.7.3. The Establishment of Objective Function

(1) Response Surface Methodology (RSM). The end effector mounting bracket is an irregular part; by the structural parameters, directly constructing accurate objective functions is difficult; therefore, Response Surface Methodology (RSM) [36] is introduced, using numerical approximation and functional fit method to construct objective functions. The basic principle of RSM is to use DOE [37] method experiments with design variables to obtain response surface model of objective functions and predict the response value of nonexperimental points. Suppose a function of design variables X and response variables F(X) can be described as(8)FX=F~X+ε=∑i=0LβiϕiX+ε,where F(X) is the actual response function of the optimization objective for the design variables, F~(X) is the approximate response function, namely, response surface model, L is the number of basis functions, ϕi(X) is a polynomial function of design variable X, and ε is the random error of the approximate value and the actual value.The approximate models of objective functions are constructed by the second-order polynomial response surface model as shown in(9)F~X=β0+∑i=1nβixi+∑i=1nβiixi2+∑i<jnβijxixj,where xi are input variables, β0 is the constant, and βi is a variable coefficient of regression model.Using the least square method solves the coefficient of (9): (10)A=β0,β1,…,βLT.Then,(11)A=XTX-1XTY,where X is the basis function matrix, which is described as (12)X=ϕ1X1ϕ2X1⋯ϕLX1ϕ1X2ϕ2X2⋯ϕLX2⋯⋯⋱⋯ϕ1XQϕ2XQ⋯ϕLXQ.(2) The Establishment of Objective Functions Based on the RSM. With minimal mass, minimizing the maximal stress, minimizing the maximal deformation, and maximizing the first natural frequency as optimization objectives, therefore objective functions are established as (13)FmX=minX∈DmX,Ff1X=maxX∈Df1X,FσX=minX∈DσmaxX,FδX=minX∈DδmaxX.According to the value range of selected design variables as well as the Box-Behnken design method generating experiment design matrix table, with a total of 41 sets, the mass (m, unit: kg); the first natural frequency (f1, unit: Hz), the maximum stress (σmax, unit: MPa), and maximum deformation (δmax, unit: mm) are calculated, respectively. The experiment design matrix and its results are obtained as shown in Tables 8 and 9, respectively.Table 8
Experiment design matrix.
Number Design variable (mm) x 2 x 3 x 4 x 6 x 8 1 120 15 16.5 140 10.5 2 120 6 16.5 140 10.5 3 80 15 16.5 140 10.5 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ 41 100 10.5 16.5 140 10.5Table 9
Experiment design results.
Number Output results m f 1 δ m a x σ m a x 1 8.0534 402.72 0.0508 5.9895 2 7.3791 399.45 0.0774 7.1030 3 8.3144 398.43 0.0489 5.9539 ⋯ ⋯ ⋯ ⋯ ⋯ 41 7.8617 399.85 0.0587 6.0313With the design variables in Table8 as inputs and the calculated results in Table 9 as outputs, the response surface models of the four objective functions are obtained as in (14)~(17), respectively, where m(X) is the mass function, f1(X) is the first natural frequency function, σmax(X) is the maximum stress function, and δmax(X) is the maximum deformation function. (14)mX=5.56-1.2×10-5x2+0.1x3+0.08x4-6.88×10-5x6+0.06x8-3.26×10-5x22-9.36×10-5x32+3.06×10-7x42+2.44×10-7x62+1.04×10-5x82+1.15×10-18x2x3+1.47×10-8x2x4+3.13×10-9x2x6+2.78×10-8x2x8-0.001x3x4+1.39×10-8x3x6+1.23×10-7x3x8+2.18×10-19x4x6+9×10-18x4x8-0.0001x6x8,(15)f1X=345.91+0.4x2+1.84x3+0.18x4+0.05x6+3.49x8-0.0004x22-0.09x32-0.01x42-0.001x62-0.03x82-0.008x2x3-0.005x2x4-0.0002x2x6+0.002x2x8+0.03x3x4+0.01x3x6-0.07x3x8-0.0013x4x6+0.007x4x8-0.002x6x8,(16)δmaxX=0.18+2.011x2-0.01x3-0.006x4+6.18×10-5x6-0.0007x8+3.08×10-7x22+0.0002x32+6.76×10-5x42+8.37×10-8x62-2.98×10-7x82-1.92×10-6x2x3-2.5×10-7x2x4+9.37×10-9x2x6-3.89×10-7x2x8+0.0002x3x4-5.58×10-6x3x6+4.4×10-5x3x8+9.26×10-7x4x6-6.34×10-6x4x8-1.39×10-8x6x8,(17)σmaxX=13.1+0.008x2-0.85x3-0.33x4+0.004x6+0.09x8+4.58×10-5x22+0.023x32+0.004x42+5.6×10-5x62+0.005x82-0.0001x2x3-0.0001x2x4-2.19x2x6-0.0009x2x8+0.016x3x4-0.0004x3x6+0.0028x3x8+8.39×10-5x4x6-0.002x4x8-0.001x6x8.In order to validate the credibility of the response surface models, reselect 10 sets’ value of design parameters in their value range by the Box-Behnken method and calculate; the calculated results were compared with response surface models predicted results as shown in Table10. Due to the fact that the mass error is 0, therefore it is not listed in the table. Table 10 shows all the errors are within 2%; the credibility of the above response surface models meets design requirements, which can be used as objective functions for the optimization design [38].Table 10
Comparison between the calculation results and predicted results of the approximate model.
Number Calculated results Response surface model predicted results Relative error (%) f 1 δ m a x σ m a x f 1 δ m a x σ m a x f 1 δ m a x σ m a x 1 398.73 0.0504 6.0405 399.05 0.0497 5.9992 0.0816 1.4429 0.6840 2 382.31 0.0552 6.2636 382.54 0.0553 6.3163 0.0591 0.1819 0.8411 3 388.72 0.0534 6.0975 389.01 0.0532 5.9920 0.0722 0.4544 1.7298 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ 10 395.56 0.0576 6.0904 395.89 0.0570 6.0237 0.0833 1.0272 1.0945


### 2.7.4. Determine the Weight Coefficients

This paper presents an analytic hierarchy process based on finite element analysis results (FEA + AHP) method to determine the weight coefficients of each optimization objective.(1) Analytic Hierarchy Process (AHP). The basic principle and steps of the AHP are as follows [39].(i) Establish Hierarchical Structure Model. Generally, a hierarchical structure model includes three levels: goal level, criteria level, and policy level, as shown in Figure 17. Hierarchy number is unrestricted; generally, the number of the elements is not more than nine in each level.Figure 17
Hierarchical structure model.(ii) Construct Pairwise Comparison Matrix. If a level has n factors, then compare their importance to a certain factor of the above level. The 1–9 scale method is used to obtain mark results as shown in Table 11. The pairwise comparison matrix A is obtained according to the mark results, with matrix A as shown in (18)A=aijn×n=a11a12⋯a1na21a22⋯a2n⋯⋯⋯⋯an1an2⋯ann.Table 11
Criteria of 1~9 scale.
Scale Scale ofi compared to j 1 Equally important 3 Weakly important 5 Strongly important 7 Very strongly important 9 Extremely important 2, 4, 6, 8 Intermediate values betweentwo adjacency values Reciprocal If marki to j is aij, then mark j to i is aij=1/aij(iii) Hierarchy Sorting and Consistency Checking. (a) Calculate the consistency index (CI): (19)CI=λmax-nn-1,where λmax is the maximum eigenvalue of the comparison matrix.(b) The random consistency index RI is shown in Table12.Table 12
Random index.
n 1 2 3 4 ⋯ 10 RI 0 0 0.58 0.90 ⋯ 1.51(c) Calculate the consistency ratio CR:(20)CR=CIRI.WhenCR<0.1, the consistency of the comparison matrix is acceptable; otherwise, the comparison matrix should be properly modified until CR<0.1.(iv) Solving the Maximum Eigenvalue λ m a x Corresponding Eigenvector W ∗ = ( ω 1 , ω 2 , … , ω n ) T of the Comparison Matrix Which Passes the Consistency Checking. And then the weight vector W is obtained by the normalization to the eigenvector W∗.(2) Determine the Weight Coefficients(i) Establishing the Hierarchical Structure Model Shown in Figure18. The goal level is structure optimization, using A to represent it. The criteria level is mainly composed of the mass, the maximum deformation, the first natural frequency, and the maximum stress, respectively, using B1, B2, B3, and B4 to represent them. The policy level is mainly composed of five design variables, respectively, using C1, C2, C3, C4, and C5 to represent them.Figure 18
Hierarchical structure model.(ii) Construct the Comparison Matrix. Based on the results of finite element analysis, the importance of each performance index is investigated. The main objective of this study is the lightweight of the end effector mounting bracket; therefore, reducing mass is the main goal. The importance of the first natural frequency, deformation, and stress is discussed. It is known from the frequency response analysis that vibration frequency corresponds to the maximum vibration velocity located in the vicinity of the first natural frequency, easily causing vibration fatigue of parts; consequently, its importance is ranked the second; as can be seen from the statics analysis, the maximum deformation is very small; however, in order to reduce the impact of the deformation on the positioning accuracy of the palletizer robot, the maximum deformation should be smaller to be better; therefore, the importance of deformation is ranked the third; the maximum stress is far less than the allowable stress; therefore, its importance comes in the last row. Accordingly, the pairwise comparison results of each factor in the criterion layer are shown in Table 13. And get the comparison matrix as shown in(21)A=1537151122132131712131.Table 13
Pairwise comparison results of each factor.
A B 1 B 2 B 3 B 4 B 1 1 5 3 7 B 2 1/5 1 1/2 2 B 3 1/3 2 1 3 B 4 1/7 1/2 1/3 1The eigenvalue of maximum of matrixA is (22)λmax=4.0192;λ m a x corresponding eigenvector is (23)W∗=0.9142,0.1912,0.3392,0.1124T.(iii) Consistency Checking Consistency index:CI=0.0064<0.1. Random consistency index:RI=0.9. Consistency ratio:CR=0.007<0.1.Accordingly, it can be determined whether comparison matrix A passes the consistency checking. Then, we carry out the eigenvector normalization process to obtain the weight vector of minimum mass, the maximum deformation, the maximum first natural frequency, and maximum stress:(24)W=0.5872,0.1228,0.2179,0.0722T.


## 2.7.1. The Choice of Design Variables

(1) Initial Choice of Design Variables. Initially, 9 structure parameters X=(x1,x2,…,x9) of the end effector mounting bracket were selected as shown in Figure 14; their names, initial values, and value range are shown in Table 5. Structural parameters are all local parameters, independent of each other, without changing the overall size and the assembly dimensions with other parts.Table 5
The initial value and range of design variables, mm.
Design variable Name Initial value Value range x 1 Small hole diameter of high upright plate 50 50~75 x 2 Large hole diameter of high upright plate 80 80~120 x 3 Arc upright plate thickness 15 6~15 x 4 Bottom plate thickness 25 8~25 x 5 Posterior slope plate thickness 15 6~15 x 6 Elliptical hole major axis of middle plate 100 100~180 x 7 Elliptical hole minor axis of middle plate 50 50~80 x 8 Middle plate thickness 15 6~15 x 9 Stiffener thickness of high upright plate 15 6~15Figure 14
Initial choice of design variables.(2) Sensitivity Analysis. In the process of optimizing design, if the design variables are too many, the key structure parameters can be neglected, resulting in the pseudo optimization results, and the calculated quantity is large and the efficiency is low. Through the sensitivity analysis, some structural design parameters, which have significantly influenced the structure response parameter, are selected as the design variable from the initial selection [30].There are many kinds of sensitivity analysis methods; among them, sensitivity analysis method based on experimental data combined with Spearman rank correlation coefficient method is a reliable and effective method [31, 32].(i) Design of Experiments (DOE) Method. Experimental data can be obtained by DOE method such as Central Composite Design (CCD), Doehlert Matrix (DM), and the Box-Behnken design (BBD) method [33]. Compared with the CCD and DM methods, Box-Behnken design method is more efficient and economical; therefore, Box-Behnken method [34] is used to estimate coefficient of the second-order polynomial response surface model. The Box-Behnken design method is an incomplete three-level part factorial experiment design method, is of nearly rotatable symmetry, is composed of multiple orthogonal cubes, and includes a central point. Schematic diagram of the three-level Box-Behnken experimental design plan is shown in Figure 15. The advantage of the Box-Behnken design is that it can avoid emergence extreme point [35].Figure 15
Schematic diagram of Box-Behnken experimental design plan.(ii) Sensitivity Analysis Based on DOE Combined with Spearman Rank Correlation Coefficient Method. Suppose sensitivity of structural response y (such as maximum stress) to structural design parameters x1,x2,x3,…,xm is examined. According to the initial value and value range of structural design parameters, we use the design of experiments (DOE) method to obtain a set of sample points of these structural design parameters x1i,x2i,x3i,…,xmi (i=1,2,…,n), and n response values of maximum stress y1,y2,y3,…,yn are obtained by n times calculating simulations. m data pairs consist of the sample values xj1,xj2,xj3,…,xjn that are calculated n times by the jth (j=1,2,3,…,m) structural design parameters, and the response values of maximum stress y1,y2,y3,…,yn by calculating simulations n times are as follows:(2)xj1y1,xj2y2,…,xjnynj=1,2,…,m.The correlation between thejth structural design parameters and the maximum stress is investigated by the Spearman rank correlation coefficient method. Respectively, the n sample values xj1,xj2,xj3,…,xjn of the jth structural design parameters and the response value y1,y2,y3,…,yn of the maximum stress are sorted in descending order, respectively, using Rji,Qi to indicate the position number of the original data xji, yi (i=1,2,…,n) in the sorted list, where Rji,Qi are called the rank of xji, yi A set of new data pairs are formed by Rji, Qi instead of xji, yi as follows: (3)Rj1Q1,Rj2Q2,…,RjnQn.Then Spearman rank correlation coefficient of the original data pairs is(4)rsj=∑i=1nRji-R-Qi-Q-∑i=1nRji-R-2∑i=1nQi-Q-2,j=1,2,3,…,m,where (5)R-=∑i=1nRjin=n+12,Q-=∑i=1nQin=n+12.Simplification of formula (1):(6)rsj=1-6∑i=1nRji-Qi2nn2-1.Rank correlation coefficientrsj, in the range of [-1,1], reflects the correlation between the maximum stress response y and the jth structural parameters. That rsj has a positive value indicates a monotonically increasing relationship between the two; that rsj=1 indicates perfect positive correlation; that rsj has a negative value indicates a monotonically decreasing relationship between the two; that rsj=-1 indicates perfect negative correlation; the greater the absolute value rsj, the greater the correlation, that is, the more sensitive the structure response y to the jth structural parameters; rsj=0 indicates that there is no relationship between two variables.By that analogy, the sensitivity of structural response parameters to all selected structural design parameters can be obtained.According to the value range of initial design variables as well as Box-Behnken design method generating experiment design matrix table, with a total of 97 sets, the mass (m, unit: kg), the first natural frequency (f1, unit: Hz), the maximum stress (σmax, unit: MPa), and the maximum deformation (δmax, unit: mm) are calculated, respectively. The experiment design matrix and corresponding results are obtained as shown in Tables 6 and 7.Table 6
Experiment design matrix, mm.
Number Design variables x 1 x 2 x 3 x 4 x 5 x 6 x 7 x 8 x 9 1 62.5 100 10.5 16.5 10.5 140 65 10.5 10.5 2 50 100 10.5 8 10.5 140 50 10.5 10.5 3 75 100 10.5 8 10.5 140 50 10.5 10.5 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ 97 62.5 100 15 16.5 15 140 80 10.5 10.5Table 7
Experiment design results.
Number Output result m f 1 δ m a x σ m a x 1 7.5701 391.18 0.0603 7.7425 2 7.0498 399.94 0.0765 9.6833 3 6.9479 399.60 0.0773 9.5524 ⋯ ⋯ ⋯ ⋯ ⋯ 97 8.0109 396.10 0.0499 7.4748Then, the Spearman rank correlation coefficient method is used to obtain the mass sensitivity, the first natural frequency sensitivity, the maximum stress sensitivity, and the maximum deformation sensitivity of the end effector mounting bracket, as shown in Figure16.Figure 16
Sensitivity analysis of structural design parameters to response parameters.
(a)
Mass sensitivity (b)
The maximum deformation sensitivity (c)
The maximum stress sensitivity (d)
The first natural frequency sensitivityHere, five structural design parameters, whose every kind of sensitivity value is greater than 0.2, are selected as design variables. The design variables are identified as(7)X=x2,x3,x4,x6,x8.


## 2.7.2. Constraint Condition

Constraint conditions are values range of every design variable.(1) Design Variable x 2 and x 6 Are the Diameter of the Hole. The greater the hole, the lighter the mass, therefore initial values of these holes diameter are set to the minimum value, their maximum values are less than boundary dimensions of structure.(2) Design Variables x 3 , x 4,andx8Are the Thickness of Plates. The smaller the thickness, the lighter the mass; therefore, the initial values of these plates’ thicknesses are set to the maximum value, and their minimum value is set to the minimum thickness value of cast aluminum alloy part [26], as shown in Table 5.


## 2.7.3. The Establishment of Objective Function

(1) Response Surface Methodology (RSM). The end effector mounting bracket is an irregular part; by the structural parameters, directly constructing accurate objective functions is difficult; therefore, Response Surface Methodology (RSM) [36] is introduced, using numerical approximation and functional fit method to construct objective functions. The basic principle of RSM is to use DOE [37] method experiments with design variables to obtain response surface model of objective functions and predict the response value of nonexperimental points. Suppose a function of design variables X and response variables F(X) can be described as(8)FX=F~X+ε=∑i=0LβiϕiX+ε,where F(X) is the actual response function of the optimization objective for the design variables, F~(X) is the approximate response function, namely, response surface model, L is the number of basis functions, ϕi(X) is a polynomial function of design variable X, and ε is the random error of the approximate value and the actual value.The approximate models of objective functions are constructed by the second-order polynomial response surface model as shown in(9)F~X=β0+∑i=1nβixi+∑i=1nβiixi2+∑i<jnβijxixj,where xi are input variables, β0 is the constant, and βi is a variable coefficient of regression model.Using the least square method solves the coefficient of (9): (10)A=β0,β1,…,βLT.Then,(11)A=XTX-1XTY,where X is the basis function matrix, which is described as (12)X=ϕ1X1ϕ2X1⋯ϕLX1ϕ1X2ϕ2X2⋯ϕLX2⋯⋯⋱⋯ϕ1XQϕ2XQ⋯ϕLXQ.(2) The Establishment of Objective Functions Based on the RSM. With minimal mass, minimizing the maximal stress, minimizing the maximal deformation, and maximizing the first natural frequency as optimization objectives, therefore objective functions are established as (13)FmX=minX∈DmX,Ff1X=maxX∈Df1X,FσX=minX∈DσmaxX,FδX=minX∈DδmaxX.According to the value range of selected design variables as well as the Box-Behnken design method generating experiment design matrix table, with a total of 41 sets, the mass (m, unit: kg); the first natural frequency (f1, unit: Hz), the maximum stress (σmax, unit: MPa), and maximum deformation (δmax, unit: mm) are calculated, respectively. The experiment design matrix and its results are obtained as shown in Tables 8 and 9, respectively.Table 8
Experiment design matrix.
Number Design variable (mm) x 2 x 3 x 4 x 6 x 8 1 120 15 16.5 140 10.5 2 120 6 16.5 140 10.5 3 80 15 16.5 140 10.5 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ 41 100 10.5 16.5 140 10.5Table 9
Experiment design results.
Number Output results m f 1 δ m a x σ m a x 1 8.0534 402.72 0.0508 5.9895 2 7.3791 399.45 0.0774 7.1030 3 8.3144 398.43 0.0489 5.9539 ⋯ ⋯ ⋯ ⋯ ⋯ 41 7.8617 399.85 0.0587 6.0313With the design variables in Table8 as inputs and the calculated results in Table 9 as outputs, the response surface models of the four objective functions are obtained as in (14)~(17), respectively, where m(X) is the mass function, f1(X) is the first natural frequency function, σmax(X) is the maximum stress function, and δmax(X) is the maximum deformation function. (14)mX=5.56-1.2×10-5x2+0.1x3+0.08x4-6.88×10-5x6+0.06x8-3.26×10-5x22-9.36×10-5x32+3.06×10-7x42+2.44×10-7x62+1.04×10-5x82+1.15×10-18x2x3+1.47×10-8x2x4+3.13×10-9x2x6+2.78×10-8x2x8-0.001x3x4+1.39×10-8x3x6+1.23×10-7x3x8+2.18×10-19x4x6+9×10-18x4x8-0.0001x6x8,(15)f1X=345.91+0.4x2+1.84x3+0.18x4+0.05x6+3.49x8-0.0004x22-0.09x32-0.01x42-0.001x62-0.03x82-0.008x2x3-0.005x2x4-0.0002x2x6+0.002x2x8+0.03x3x4+0.01x3x6-0.07x3x8-0.0013x4x6+0.007x4x8-0.002x6x8,(16)δmaxX=0.18+2.011x2-0.01x3-0.006x4+6.18×10-5x6-0.0007x8+3.08×10-7x22+0.0002x32+6.76×10-5x42+8.37×10-8x62-2.98×10-7x82-1.92×10-6x2x3-2.5×10-7x2x4+9.37×10-9x2x6-3.89×10-7x2x8+0.0002x3x4-5.58×10-6x3x6+4.4×10-5x3x8+9.26×10-7x4x6-6.34×10-6x4x8-1.39×10-8x6x8,(17)σmaxX=13.1+0.008x2-0.85x3-0.33x4+0.004x6+0.09x8+4.58×10-5x22+0.023x32+0.004x42+5.6×10-5x62+0.005x82-0.0001x2x3-0.0001x2x4-2.19x2x6-0.0009x2x8+0.016x3x4-0.0004x3x6+0.0028x3x8+8.39×10-5x4x6-0.002x4x8-0.001x6x8.In order to validate the credibility of the response surface models, reselect 10 sets’ value of design parameters in their value range by the Box-Behnken method and calculate; the calculated results were compared with response surface models predicted results as shown in Table10. Due to the fact that the mass error is 0, therefore it is not listed in the table. Table 10 shows all the errors are within 2%; the credibility of the above response surface models meets design requirements, which can be used as objective functions for the optimization design [38].Table 10
Comparison between the calculation results and predicted results of the approximate model.
Number Calculated results Response surface model predicted results Relative error (%) f 1 δ m a x σ m a x f 1 δ m a x σ m a x f 1 δ m a x σ m a x 1 398.73 0.0504 6.0405 399.05 0.0497 5.9992 0.0816 1.4429 0.6840 2 382.31 0.0552 6.2636 382.54 0.0553 6.3163 0.0591 0.1819 0.8411 3 388.72 0.0534 6.0975 389.01 0.0532 5.9920 0.0722 0.4544 1.7298 ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ ⋯ 10 395.56 0.0576 6.0904 395.89 0.0570 6.0237 0.0833 1.0272 1.0945


## 2.7.4. Determine the Weight Coefficients

This paper presents an analytic hierarchy process based on finite element analysis results (FEA + AHP) method to determine the weight coefficients of each optimization objective.(1) Analytic Hierarchy Process (AHP). The basic principle and steps of the AHP are as follows [39].(i) Establish Hierarchical Structure Model. Generally, a hierarchical structure model includes three levels: goal level, criteria level, and policy level, as shown in Figure 17. Hierarchy number is unrestricted; generally, the number of the elements is not more than nine in each level.Figure 17
Hierarchical structure model.(ii) Construct Pairwise Comparison Matrix. If a level has n factors, then compare their importance to a certain factor of the above level. The 1–9 scale method is used to obtain mark results as shown in Table 11. The pairwise comparison matrix A is obtained according to the mark results, with matrix A as shown in (18)A=aijn×n=a11a12⋯a1na21a22⋯a2n⋯⋯⋯⋯an1an2⋯ann.Table 11
Criteria of 1~9 scale.
Scale Scale ofi compared to j 1 Equally important 3 Weakly important 5 Strongly important 7 Very strongly important 9 Extremely important 2, 4, 6, 8 Intermediate values betweentwo adjacency values Reciprocal If marki to j is aij, then mark j to i is aij=1/aij(iii) Hierarchy Sorting and Consistency Checking. (a) Calculate the consistency index (CI): (19)CI=λmax-nn-1,where λmax is the maximum eigenvalue of the comparison matrix.(b) The random consistency index RI is shown in Table12.Table 12
Random index.
n 1 2 3 4 ⋯ 10 RI 0 0 0.58 0.90 ⋯ 1.51(c) Calculate the consistency ratio CR:(20)CR=CIRI.WhenCR<0.1, the consistency of the comparison matrix is acceptable; otherwise, the comparison matrix should be properly modified until CR<0.1.(iv) Solving the Maximum Eigenvalue λ m a x Corresponding Eigenvector W ∗ = ( ω 1 , ω 2 , … , ω n ) T of the Comparison Matrix Which Passes the Consistency Checking. And then the weight vector W is obtained by the normalization to the eigenvector W∗.(2) Determine the Weight Coefficients(i) Establishing the Hierarchical Structure Model Shown in Figure18. The goal level is structure optimization, using A to represent it. The criteria level is mainly composed of the mass, the maximum deformation, the first natural frequency, and the maximum stress, respectively, using B1, B2, B3, and B4 to represent them. The policy level is mainly composed of five design variables, respectively, using C1, C2, C3, C4, and C5 to represent them.Figure 18
Hierarchical structure model.(ii) Construct the Comparison Matrix. Based on the results of finite element analysis, the importance of each performance index is investigated. The main objective of this study is the lightweight of the end effector mounting bracket; therefore, reducing mass is the main goal. The importance of the first natural frequency, deformation, and stress is discussed. It is known from the frequency response analysis that vibration frequency corresponds to the maximum vibration velocity located in the vicinity of the first natural frequency, easily causing vibration fatigue of parts; consequently, its importance is ranked the second; as can be seen from the statics analysis, the maximum deformation is very small; however, in order to reduce the impact of the deformation on the positioning accuracy of the palletizer robot, the maximum deformation should be smaller to be better; therefore, the importance of deformation is ranked the third; the maximum stress is far less than the allowable stress; therefore, its importance comes in the last row. Accordingly, the pairwise comparison results of each factor in the criterion layer are shown in Table 13. And get the comparison matrix as shown in(21)A=1537151122132131712131.Table 13
Pairwise comparison results of each factor.
A B 1 B 2 B 3 B 4 B 1 1 5 3 7 B 2 1/5 1 1/2 2 B 3 1/3 2 1 3 B 4 1/7 1/2 1/3 1The eigenvalue of maximum of matrixA is (22)λmax=4.0192;λ m a x corresponding eigenvector is (23)W∗=0.9142,0.1912,0.3392,0.1124T.(iii) Consistency Checking Consistency index:CI=0.0064<0.1. Random consistency index:RI=0.9. Consistency ratio:CR=0.007<0.1.Accordingly, it can be determined whether comparison matrix A passes the consistency checking. Then, we carry out the eigenvector normalization process to obtain the weight vector of minimum mass, the maximum deformation, the maximum first natural frequency, and maximum stress:(24)W=0.5872,0.1228,0.2179,0.0722T.


## 2.8. Solving Multicriteria Optimization Problem

### 2.8.1. Elitist Nondominated Sorting Genetic Algorithm (NSGA-II)

NSGA-II algorithm is an improved algorithm of NSGA, which was proposed by Deb et al. [40]. NSGA-II algorithm uses fast nondominated sorting procedure, elite strategy, and nonparameter ecological niche operator, overcoming the shortcomings of the traditional NSGA such as high computational complexity, nonelite strategy, and the need for specially designated shared radius. NSGA-II algorithm has been widely used due to its good performance in exploring; Pareto set, which is obtained by NSGA-II algorithm, has good accuracy and dispersity. The NSGA-II procedure is shown in Figure 19 [7, 40].Figure 19
The main flow chart of NSGA-II.Step 1.
Initially, a random parent populationP0 is created. The population is sorted based on the nondomination. Each solution is assigned a fitness (or rank) to its nondomination level; binary tournament selection, recombination, and mutation operators are used to create a child population Q0 of size N and set t=0.Step 2.
A combined populationRt=Pt∪Qt is formed. Here, Pt is the parent population in tth iteration; Qt is the child population in tth iteration; t denotes the iteration or generation number; the population Rt is of size 2N. Then, Rt is sorted according to nondomination. Since all previous and current population members are included in Rt, elitism is ensured. And then nondominated fronts F1,F2,…,Fi are obtained.Step 3.
SortFi using the crowded-comparison operator <n in descending order and choose the best solutions to form a new population Pt+1.Step 4.
The new populationPt+1 of size N is used for selection, crossover, and mutation to create a new population Qt+1 of size N.Step 5.
If the termination conditions are satisfied, it is over; otherwise,t=t+1; go to Step 2.Parameter configuration of NSGA-II algorithm: population size is 40; number of generations is 200; crossover rate is 0.9; cross distribution index is 10; mutation distribution index is 20.Justification of NSGA-II algorithm parameters selection: firstly, population size directly affects astringency or computational efficiency. If population size is too small, optimization results are likely to converge to local optimal solution. If population size is too large, the calculations will take a long time. Therefore, population size usually is selected in the range 20~200. Secondly, the meaning of number of generations is that NSGA-II algorithm will be terminated when the number of evolutionary generations reaches this set value. Generally, the number of generations is selected in the range 100~500 [41]. In this work, the above two parameters are determined by some tests for permutation and combination of different values of population size and number of generations in their respective ranges. The comparison test results indicate that there is only small difference in the decimal part among the group results; however, precision of results is improved with the increase of population size and number of generations, and the calculation time has become longer. Population size = 40 and number of generations = 200 are determined by taking into account the calculation precision and time. Thirdly, crossover rate controls frequentness of crossover operation. If crossover rate is too high, generation gap is likely to happen. However, if crossover rate is too low, many individuals will be prone to duplicating in the next generation directly, leading to stagnation during searching. Therefore, crossover rate is usually selected in the range 0.4~0.99 [41]. In some literatures, crossover rate is 0.9 [38, 42]. Finally, cross distribution index and mutation distribution index both are nonnegative number and user-defined rule, and the higher their values, the smaller the probability of offspring individuals far away from parent individuals, and vice versa. Therefore, the selected cross distribution index and mutation distribution index, respectively, are 10 and 20; they are smaller, in order to enhance searching ability of algorithm [38, 43].


## 2.8.1. Elitist Nondominated Sorting Genetic Algorithm (NSGA-II)

NSGA-II algorithm is an improved algorithm of NSGA, which was proposed by Deb et al. [40]. NSGA-II algorithm uses fast nondominated sorting procedure, elite strategy, and nonparameter ecological niche operator, overcoming the shortcomings of the traditional NSGA such as high computational complexity, nonelite strategy, and the need for specially designated shared radius. NSGA-II algorithm has been widely used due to its good performance in exploring; Pareto set, which is obtained by NSGA-II algorithm, has good accuracy and dispersity. The NSGA-II procedure is shown in Figure 19 [7, 40].Figure 19
The main flow chart of NSGA-II.Step 1.
Initially, a random parent populationP0 is created. The population is sorted based on the nondomination. Each solution is assigned a fitness (or rank) to its nondomination level; binary tournament selection, recombination, and mutation operators are used to create a child population Q0 of size N and set t=0.Step 2.
A combined populationRt=Pt∪Qt is formed. Here, Pt is the parent population in tth iteration; Qt is the child population in tth iteration; t denotes the iteration or generation number; the population Rt is of size 2N. Then, Rt is sorted according to nondomination. Since all previous and current population members are included in Rt, elitism is ensured. And then nondominated fronts F1,F2,…,Fi are obtained.Step 3.
SortFi using the crowded-comparison operator <n in descending order and choose the best solutions to form a new population Pt+1.Step 4.
The new populationPt+1 of size N is used for selection, crossover, and mutation to create a new population Qt+1 of size N.Step 5.
If the termination conditions are satisfied, it is over; otherwise,t=t+1; go to Step 2.Parameter configuration of NSGA-II algorithm: population size is 40; number of generations is 200; crossover rate is 0.9; cross distribution index is 10; mutation distribution index is 20.Justification of NSGA-II algorithm parameters selection: firstly, population size directly affects astringency or computational efficiency. If population size is too small, optimization results are likely to converge to local optimal solution. If population size is too large, the calculations will take a long time. Therefore, population size usually is selected in the range 20~200. Secondly, the meaning of number of generations is that NSGA-II algorithm will be terminated when the number of evolutionary generations reaches this set value. Generally, the number of generations is selected in the range 100~500 [41]. In this work, the above two parameters are determined by some tests for permutation and combination of different values of population size and number of generations in their respective ranges. The comparison test results indicate that there is only small difference in the decimal part among the group results; however, precision of results is improved with the increase of population size and number of generations, and the calculation time has become longer. Population size = 40 and number of generations = 200 are determined by taking into account the calculation precision and time. Thirdly, crossover rate controls frequentness of crossover operation. If crossover rate is too high, generation gap is likely to happen. However, if crossover rate is too low, many individuals will be prone to duplicating in the next generation directly, leading to stagnation during searching. Therefore, crossover rate is usually selected in the range 0.4~0.99 [41]. In some literatures, crossover rate is 0.9 [38, 42]. Finally, cross distribution index and mutation distribution index both are nonnegative number and user-defined rule, and the higher their values, the smaller the probability of offspring individuals far away from parent individuals, and vice versa. Therefore, the selected cross distribution index and mutation distribution index, respectively, are 10 and 20; they are smaller, in order to enhance searching ability of algorithm [38, 43].


## 3. Calculation Results and Analyses

Recommended values of structure parameters are obtained by the NSGA-II algorithm, and then they are adjusted in a small range with comprehensive consideration of casting process parameters, mounting holes on the influence of the actual structure, and many simulation tests with changing parameters; the optimization results are shown in Table14, using the optimized structure parameters to regenerate 3D model. Structure response parameters of the optimized model are obtained by static analysis and modal analysis and compared with structure response parameters of the initial model as shown in Table 15. The comparison results indicate that, after the structure optimization, the mass is reduced by 17.8%; the maximum deformation is increased by about 0.012 mm, and the maximum stress value is increased by about 1.7 MPa, but still far less than the allowable value; the first natural frequency is increased by about 9.5 Hz. And then the frequency response analysis for the optimized model is conducted again, and the analysis result indicates that the maximum vibration velocity still occurred in the same region as node A; however, the #1 peak vibration velocity at 430 Hz almost completely disappeared, moreover, the #2 peak vibration velocity is reduced by about 75%, significantly reducing the influence of vibration on the reliability and motion accuracy of the palletizing robot. The stress nephogram, the deformation nephogram, the first modal shape of optimized model, and comparison of frequency response curve for the optimized model and the initial model are shown in Figures 20–24, respectively.Table 14
Optimization results of structural parameters, mm.
Structural parameter x 2 x 3 x 4 x 6 x 8 Initial value 80 15 25 100 15 Recommended value 114.7 10.87 15.08 100.08 14.99 Optimized value 100 10 16 100 15Table 15
Optimization results of objective parameters.
Objective parameters m f 1 δ m a x σ m a x Initial value 9.790 409.6 0.042 6.806 Optimized value 8.048 419.1 0.054 8.473 Variation −1.742 9.5 0.012 1.667Figure 20
Stress nephogram of optimized model.Figure 21
Deformation nephogram of optimized model.Figure 22
The first modal shape of optimized model.Figure 23
Frequency response curve of optimized model.Figure 24
Comparison of frequency response curve for the optimized model and the initial model.


## 4. Conclusion

(1) In this paper, the multicriteria structure optimization problem of the end effector mounting bracket of MD-1200YJ type palletizing robot was researched. The finite element model of the part is established, and the accuracy of the finite element model was verified by modal test.(2) Through modal analysis, vibration response test, frequency response analysis, and the static analysis, considering inertia, the mass is minimized, the maximal stress is minimized, the maximal deformation is minimized, and the first natural frequency is maximized as the optimization objectives are determined. Five design variables were selected by sensitivity analysis based on the DOE combined with the Spearman rank correlation coefficient method, as well as their value range as the constraint conditions; the Box-Behnken method was used to construct the RSM model of objective functions and verify the credibility of the RSM models meeting the requirements, putting forward an analytic hierarchy process based on the finite element analysis (FEA + AHP) method, to improve the objectivity of comparison matrix to a certain extent, to determine the weight coefficients of objective functions; the multicriteria optimization design model was established and the NSGA- II algorithm is used for solving to obtain optimization results.(3) Compared with the initial model, mass of the optimized model is reduced by 17.8% in the case of the maximum stress and deformation within the allowable range; and the first natural frequency of the optimized model is increased and vibration response characteristics of the entire structure are improved significantly. The validity of the optimization design method is verified.



---

*Source: 6049635-2018-04-22.xml*