Nonlinear friction dynamic modeling and performance analysis of flexible parallel robot

The article presents a friction dynamic modeling method for a flexible parallel robot considering the characteristics of nonlinear friction. An approximate friction model which is proposed by Kostic et al. is applied to establish the dynamic model with Lagrange method. Parameters identification is completed by least square method, and the tracking accuracy and the stability of the robot are systematically analyzed before and after dynamic compensation at different speeds. Its position error of the robot after compensation is only 0.98 mm at low speed. The accuracy is improved 10 times than that before compensation. In addition, the peak velocity errors are 3.97 mm·s−1 and 1.49 mm·s−1 at high and low speed, respectively, which are reduced 2.5 times than that before compensation. The experimental data also indicate that velocity tracking curve has no obvious peak error compared with the common method based on Coulomb and viscous friction model. The curve is much smoother with proposed model, and the motion stability of robot at high speed has been greatly improved. The proposed method just needs the robot to collect some positions before path tracing, and the parameters identification of dynamic model can be completed quickly. The compensation effect is much more better than common method. So the proposed method can be extended to complete the dynamic identification for complex robot with more joints. It is helpful to further improve the stability and the accuracy at high speed.


Introduction
The development of robot technology has great significance in improving the competitiveness of advanced production and promoting the national economy. The serial robot has the characteristics of simple structure and easy operation. It was first used in industrial production. [1][2][3][4] However, it has some problems of error accumulation effect and low load because of single chain structure. The parallel mechanism has the advantages of strong load capacity, large stiffness, and no accumulation of errors. 5,6 But it is difficult to build dynamic model because of multiple chains and a large number of joints. 7,8 At present, kinematic calibration and error compensation are usually used to improve the accuracy of robot. For improving the stability of the mechanism, the research on dynamic parameter identification has important theoretical research engineering. [9][10][11][12] The dynamic model of the parallel robot can describe the relationship between the robot's motion and the torque of each joint. The stability of the robot's actual trajectory tracking can be reflected by the friction force curve. It has proved that, the nonlinear friction of the driving joints makes an increasingly prominent impact on the accuracy, motion stability, and rapidity of the robot. 13 The dynamic modelling methods of parallel robot can be divided into Newton-Euler method and Lagrange method. Compared with Newton-Euler method, the dynamic model which is obtained by Lagrange method is more concise. Therefore, in the field of parallel robot design, Lagrange method has been widely used. For the high-speed parallel mechanism, the inertia load of each linkage is directly proportional to rotating speed of the original component. With the increase of velocity, the effect of joint friction becomes more obvious, and it often results in different degrees of shaking on the trajectory tracking. 14 Therefore, the dynamic identification and compensation to improve the stability of robot system have gradually become a focus of attention. The friction is the main factor that affects the stability of the robot system. It will cause large steady-state error of the system and a dead zone near zero speed. 15 Therefore, scholars have carried out research on the speed control strategy. 16,17 At present, the main velocity planning algorithm includes ladder velocity control algorithm (LVCA) and Stype velocity control algorithm (SVCA). The latter is named for its S-shaped speed curve, which is widely used in high-speed and high-precision machining systems. 15,18 The SVCA is suitable for 3D trajectory tracking of robot with high stability and accuracy, but it owns many parameters which result the algorithm become more complex. 19,20 The LVCA is relatively simple and easy to be implemented. The accuracy is poorer than SVCA when robot moves on 3D trajectory path. So the LVCA is widely applied in the control of planar industrial robots. The researched robot in the article is mainly used for precise positioning and arbitrary trajectory tracking in the work plane. It is widely used in the industrial field because of its low cost, simple structure, and easy control.
Many scholars have discussed and summarized the friction characteristics and models. Berger 21 introduced the friction models and typical dynamic response in different fields. Liu et al. 22 introduced the friction phenomenon and common friction model in mechanical system, which provided reference for the selection of friction model in the research, and the results have proved that the compensation effect only based on Coulomb friction, Stribeck friction, viscous friction, and static friction model is not good for robot. The main reason is that the above models are difficult to accurately and comprehensively reflect the nonlinear characteristics of friction force under dynamic motion at high speed. Some researchers focus on the influence of Coulomb friction and viscous friction on the motion trajectory. They propose a convex relaxation method to complete the optimal trajectory and feed-forward control, which is applied to solve the nonconvex problem of the manipulator. 13,23 But this method does not consider the effect of Stribeck friction. Some studies only consider the influence of the viscous friction effect on the motion trajectory. 24 However, many researchers just consider the influence of Coulomb friction and viscous friction on the trajectory of the robot, which ignores the actual nonlinear distribution of the manipulator during acceleration and deceleration process. 10,11 Kostic et al. 25 proposed an approximate nonlinear friction model, which reflects the coupling characteristics of the Coulomb friction, viscous friction, and Stribeck friction. The exponential function is used to reflect the nonlinear friction characteristic of driving joints during acceleration and deceleration, which is verified in the RRR manipulator. The dynamic model of the robot considering nonlinear friction is established with Lagrange method. An approximate exponential friction model is introduced to dynamic model of robot, and the parameter identification of dynamic is completed by least square method (LSM). On this basis, the robot's dynamic performance has been analyzed deeply before and after compensation. The feasibility and effectiveness of proposed method are verified by experiments.
In this article, a complete and systematic identification process is presented for an industrial parallel robot using an approximate friction model and LSM. The proposed method just needs to collect some positions data. But these sampling position points should be distributed in the whole working space as far as possible, and the sampling accuracy must be ensured when robot is moving on the preset path. So the article recommends to apply single point sampling method. The proposed method can be used to finish dynamic identification work for other robots with more complex structure. The rest of this article is organized as follows. The second section presents the nonlinear friction dynamic modeling method, which includes the robot's structure, Lagrange modeling principle and the adopted exponential friction model. The third section includes experimental set, sampling method of position, stability analysis, and validation of proposed method. The last section presents the conclusions.

Nonlinear friction dynamic modeling
The structure of parallel robot The researched robot in this article is produced by Shenzhen GuGao Tech Company in China, which is mainly used for precise positioning and trajectory tracking on the working plane. It owns some advantages of low cost, simple structure, and easy control. Therefore, the robot has been widely used in industrial field. 14,17 A large number of scholars have made much research on the kinematic modeling, structure error calibration, and dynamic modeling about plane parallel robot. But there are few research on the dynamic parameter identification and compensation effect under considering nonlinear friction.
Fu et al. 13 has demonstrated that friction force of the drive joints is an important factor that affects the robot accuracy and movement stability on continuous path. In order to improve the robot's tracing accuracy and stability at high speed, the parameters identification method of dynamics is researched deeply and systematically considering nonlinear friction factor. The robot's structure is shown in Figure 1. F i ðx Di ; y Di Þ (i ¼ 1, 2, 3) means the position of the driving motor which is fixed on the platform in Figure 2, and G i ðx Bi ; y Bi Þ (i ¼ 1, 2, 3) is the passive joint's coordinate. L Di and L Bi are respectively the length of driving rod and the passive rod. Oðx o ; y o Þ represents the robot's coordinate. q i ði ¼ 1; 2; 3Þ and q Bi ði ¼ 1; 2; 3Þ are called the driving angle and the passive angle, respectively.
So the coordinate of the robot can be written as with Therefore, the following conclusion is deduced as Then the equation (3) can be modified as where C i ði ¼ 1; 2; 3Þ are defined as Therefore, the coordinate of robot can be written as and the speed of robot is also expressed as So the velocity _ q can be written as where the matrix J is defined as  Figure 2. The robot's coordinate system.
All elements of matrix J are defined as following

Lagrange dynamic modeling
The Lagrange function L can be written as 7 where T represents kinetic energy of system and U denotes the potential energy of system. Because this robot works on the xoy plane, the potential energy U of system will be ignored. The structure of single chain is given in Figure 3. D i ðx cdi ; y cdi Þ and B i ðx cbi ; y cbi Þ are respectively gravity center of the driving rod and the passive rod, respectively. So, the kinetic energy of the driving rod and the passive rod can be written as where _ q represents the driving joint's velocity. I D and I B are called moment of inertia. m D and m B are the quality of the driving rod and the passive rod, respectively, and x cdi , y cdi , x cbi , and y cbi are defined as with where r d and r b represent the length of driving rod and the length of passive rod. Therefore, the Lagrange function of the robot is expressed as In which these parameters e i , i , and k i are defined as Then the Lagrange equation can be expressed as where t means the torque of the joints, it is defined as So the dynamic model can be written as The above equation (26) can be simplified as Figure 3. The structure of single chain.
In which M and C are defined as following The constraint condition of the researched robot is written as HðqÞ ¼ H f 1 ðqÞ f 2 ðqÞ ::::: f m ðqÞ where the matrix HðqÞ is defined as HðqÞ ¼ x D1 þ l cos q 2 þ l cos q b1 À x D2 À l cos q 2 À l cos q b2 y D1 þ l sin q 2 þ l sin q b1 À y D2 À l sin q 2 À l sin q b2 x D1 þ l cos q 2 þ l cos q b1 À x D3 À l cos q 3 À l cos q b3 y D1 þ l sin q 2 þ l sin q b1 À y D3 À l sin q 3 À l sin q b3 Then the velocity constraint of any point in the workspace can be expressed as In this formula, SðqÞ is called the velocity constraint matrix, and the acceleration constraint condition also can be obtained as 7 Therefore, the dynamic model of system including the velocity constraint can be obtained as where S T is the velocity constraint matrix, and h is defined as Àl sin q b1 l sin q b2 0 l cos q 1 Àl cos q 2 0 l cos q b1 Àl cos q b2 0 Àl sin q 1 0 l sin q 3 Àl sin q b1 0 l sin q b3 l cos q 1 0 Àl cos q 3 l cos q b1 0 Àl cos q b3 Then the dynamic model can be further modified as where M s , C s , and t s are defined as The equation (36) can be also written as where ðJ T Þ þ is the inverse matrix of J T .

Nonlinear friction modeling
Fu et al. and Yang et al. 13,26 have proved that the friction force of the driving joints is the main factor which affects the robot's accuracy and stability under continuous path tracing. Therefore, it is very important to establish an accurate nonlinear friction dynamic model. An approximate frictional force model which synthetically considers common friction characteristics such as Coulomb and Stribeck friction has been used to research on dynamic modeling. It is proposed by Hensens and Kosticwhich is defined as 12 where _ q j i is the i joint's angular velocity at the j position of path, and k vi represents viscous friction coefficient. m 1i and m 2i are the approximation coefficient for Coulomb and Stribeck friction. f 1i and f 2i represent the Coulomb and the Stribeck friction which are defined as follows Therefore, the dynamic model considering active joints' friction can be written as Because the input torque of control system is digital, the ratio between the standard unit of torque and the quantitative torque is defined as k f . The nonlinear friction dynamic model can be expressed as where E, M s , and C s are defined as follows Therefore, the optimal function can be written as where t j qi and E j i are the i joint's driving torque and dynamic torque at the j position, respectively. Sixteen parameters such as k vi , f 1i , f 2i , m 1i , m 2i , and k f (i ¼ 1, 2, 3) will be identified which are dimensionless parameters. We need to calculate the minimum value of optimization function f, and t j qi , E j i , and _ q j i can be obtained from robot system. We will apply the function of LSM which is provided by MATLAB2016a to complete the optimization calculation based on experimental data. The LSM has some advantages of high efficiency and mature algorithm which are widely used in solving calculation of extremum problem. A set of parameters will be obtained which can make the objective complete the above sampling work under LVCA at a speed of 300 mmÁs À1 .

Identification of dynamic parameters
The single-point sampling style is selected to carry out the experiment. The robot stays 2 s at each point to ensure the reliability and accuracy of data collection. All the sampling positions are shown in Figure 4. The 10 green points are distributed on boundary and 11 red points are close to the center area.
The robot's structure parameters are given in Table 1, and the result of identification is shown in Table 2.
The curve of location error curve data can be obtained before and after compensation which is given in Figure 5. It indicates that the position error of the working space boundary (Number 0 to 9) is much larger than the center area. Before compensation, the maximum error is up to 7.09 mm at the boundary area, and it is reduced to 1.88 m after compensation. In addition, the location error in the center area is just ranged from À0.18 mm to 0.19 mm after compensation. Therefore, the nonlinear friction compensation can greatly improve location accuracy in the whole working space. In addition, Figure 5 also indicates that the robot owns better working accuracy in the center area before and after compensation.
The above data in Figure 5 are further solved and its corresponding histogram including maximum error, average error, and minimum error can be obtained shown in Figure 6.  It is not difficult to find that the maximum error of the robot appears in the boundary area. The maximum error before compensation is 7.09 mm, which is reduced to 1.88 mm after compensation. The minimum errors are distributed in the center area which are 0.45 and 0.13 mm, respectively. The average error indicates that the location error is 3.88 mm before compensation. After dynamic compensation, the error is only 0.55 mm which is reduced by an order of magnitude. So the robot's control accuracy is greatly improved on the whole path.

Experiment setup
In order to verify the effectiveness of proposed method, the robot's dynamic performance will be tested on the continuous motion. In order to ensure the safety of the robot, it is necessary to define the maximum speed and acceleration of the system when we begin to plan the path. Cai 7 provides the constraint conditions of the maximum velocity V max and the maximum acceleration a max in workspace W. They are defined as where r max ðJ þ Þ is the maximum singular value of J þ matrix which can be quickly calculated based on dynamic equation by MATLAB, and r max ðM s À1 J T Þ represents the largest singular value of matrix M s À1 J T . jjt q jj max is the maximum driving torque which can be obtained from motor product manual. n means the reduction ratio, and ! max is called the maximum angular velocity of motor. n and ! max are considered to be 100 and 3500 rÁmin À1 , respectively. The motion path is set as a line shown in Figure 7, and E 1 and E 2 are called the start point and the end point, respectively. The postion coordinates are shown in Table 3. Figure 8 is the control curve of LVCA. Table 4 has given the maximum velocity and the maximum acceleration values. In order to ensure the motion security of the robot, actual control parameters are set as 0.40 and 7.00 m Á s À2 . The experimental steps are concluded as follows; Step (1): E 1 ðx q ; y q Þ and E 2 ðx e ; y e Þ are coordinate of the start point and the end point, respectively. represents forward angle which can be calculated by ¼ arctanðy e À y q ; x e À x q Þ.
Step (2): v max and a max are given in Table 4, and robot is controlled under LVCA.
Step (3): sðtÞ is the path function. So the theoretical position at any time can be calculated as follows Step (4): The theoretical velocity and acceleration can be easily obtained at X and Y directions, respectively v x ðtÞ ¼ vðtÞcos ; v y ðtÞ ¼ vðtÞsin a x ðtÞ ¼ aðtÞcos ; a y ðtÞ ¼ aðtÞsin ð49Þ Step (5): The sampling time T c is set 0.001 s, and i is the number of sampling points; the motion time T is quickly calculated with the control parameters. Then the number of sampling points is easy to be com- Step (6): The actual speed of the robot can be obtained by sensors. Then the velocity tracking error at any time is defined as Therefore, the angle velocity curve of each driving joint can be obtained based on the above steps. Figure 9 indicates that when the angular velocity is low, the difference before and after compensation is not obvious about the driving joint 3. It can be seen that the velocity of the joint 3 during uniform stage emerges suddenly changes compared with the other two joints, which is up to 1.81 radÁs À1 without friction compensation shown in Figure 9(a). It also proves that the dynamic performance of the motor is much poorer and friction force seriously affects the stability of the robot's motion. Figure 9(b) indicates that the stability of all driving joints is greatly improved after compensation. The maximum sudden change of velocity is just 0.58 radÁs À1 which is three times smaller than before compensation. Each joint's velocity curve on the acceleration stage, the uniform stage, and the deceleration stage are given in Figures 10 to 12, respectively.  Figures 10 to 12 clearly present that motors can't work at a constant value in the stage of uniform velocity. The main reason is that nonlinear friction force leads to speed loss in different degrees, and the friction force is closely related to the motor speed. In order to research the coupling mechanism, the nonlinear friction force of each driving joint has been measured. According to the feedback information from speed sensors of driving joints, the friction force curve can be quickly drawn based on the adopted model in this article. Figure 11 indicates that the maximum friction force is distributed in the range of 2.0-7.4 N when the robot is working at uniform stage.  The force variation of driving joint 2 is 0.7 N in the stage of uniform velocity, which is the biggest variation among the three joints. It also demonstrates the stability of driving joint 2 is relatively poorer than other two joints. It shows that nonlinear friction compensation can't eliminate the influence on the stability of the robot in theory. But accurate model is more helpful to improve the stability by compensation. The maximum value of friction force variation is just 0.7 N at high speed, which indicates the three driving motors work on better condition after compensation. The two exponential functions were applied to the adopted model of friction force to approximate Coulomb and Stribeck friction. The friction force curve of joint 2 shows the characteristics of Coulomb and Stribeck friction are more obvious than viscous friction. It makes the curve appearing exponential distribution in the process of acceleration and deceleration. It can be seen that the friction curve can reflect the characteristic of the speed and the motion stability of the system in different stages. Therefore, in the actual work, the stability of the robot can be quickly analyzed by the velocity tracking performance. Experiment results are given in Figure 13 and Table 5. Figure 13 indicates that the friction force curves of the three driving joints are smooth and the force amplitude is just distributed in the range of 2.2N-7.5N after compensation. It proves the motion stability of robot is good. Table 5 demonstrates that the velocity tracing error becomes greater with the increase of motion speed. When robots work at 400 mm s À1 , the maximum and  minimum velocity peak errors are 9.80 mm s À1 and À9.99 mm s À1 before compensation, respectively. The data indicate that the stability of the robot is much poorer during acceleration and deceleration process, and its average speed error is up to 7.76 mm s À1 . It proves that the tracing accuracy is much poorer than motion on low speed which is consistent with the actual situation. After compensation, the performance of speed tracking is improved obviously. The maximum and minimum velocity peak errors are reduced to 3.97 mm s À1 and À3.86 mm s À1 , respectively. It shows that the accuracy and the stability of the motion are greatly improved on the whole path. In addition, when the robot's speed is closer to the maximum speed limit of 0.42 m s À1 , the stability and accuracy become much poorer before and after compensation. The experimental data have proved that the peak speed error and the average speed error reduce obviously with the speed decreases. The performance of velocity tracing is the best when robot works at 100 mm s À1 based on the four-group experimental data. The maximum speed peak error and the average speed error are 3.98 and 1.94 mm s À1 before compensation, which become to be 1.49 and 1.32 mm s À1 after compensation. The above results also indicate that when the speed is lower, the dynamic compensation effect is not obvious before and after compensation. Its average error difference is only 0.17 mm s À1 before and after compensation. In order to further improve the position accuracy of the robot, kinematic compensation is needed to be completed to further improve accuracy.

Compensation validation
In order to verify the effectiveness of dynamic compensation, we make an experiment to analyze the tracing  performance compared with the common friction model. Coulomb and viscous friction force model can be expressed as where k v , v, and f c are viscosity coefficient, velocity of motor, and Coulomb friction force, respectively. sgnðvÞ is defined as 22 Then the velocity tracing error of the robot is tested before and after compensation at different speeds. Figure 14 indicates that the velocity curve has many peaks and tracing error become larger with the speed increase before compensation. The motion stability of the robot is much poorer at high and low speed because of nonlinear friction force.
The smoothness of velocity curve and the tracking accuracy are lower than proposed method in this article at high speed based on Coulomb and viscous model. It is difficult to accurately reflect the nonlinear friction characteristics at high speed based on the linear friction model. The adopted approximate model includes viscous friction and uses two exponential functions to approximate the Coulomb friction and the Stribeck friction. It can truly express the robot's dynamic characteristics of nonlinear friction at high speed. Therefore, the velocity tracing accuracy and stability are much better than the common method.

Conclusions
In this article, an approximate exponential friction model is applied to establish the nonlinear friction dynamic model for a planar parallel robot. The work of dynamic identification can be completed with LSM just by measuring some position data. The maximum position error is 1.88 mm, which is seven times lower than that before compensation. The average error after compensation is only 0.55 mm, which is reduced an order of magnitude than that before compensation. So friction compensation can improve the location accuracy. Four-group experiments on continuous trajectory tracking have proved that the speed tracking accuracy and stability of the robot at high speed are worse than that at low speed. Before compensation, the velocity curve obviously has many peaks at high speed. Its maximum speed peak error is up to 9.8 mm s À1 , and the average error is 7.76 mm s À1 . It indicates the motion stability is much poorer in the process of acceleration and deceleration transition. After compensation, the maximum error and average error are reduced to 3.97 and 3.43 mm s À1 , and the velocity curve is more smoother, which shows that robot's stability and the tracing performance are greatly improved, and the friction force curve also confirms this conclusion. The average speed error is 1.94 and 1.32 mm s À1 , respectively, before and after compensation when robot works at a low speed of 100 mm s À1 . The effect of dynamic compensation is not obvious. Although the stability and velocity tracking accuracy are obviously better than that at high speed, the rapidity of competing tracking is also reduced. The velocity curve of the robot is obviously smoother than that based on Coulomb and viscous model. The results also prove that the adopted friction model in the article can be more realistic to reflect the characteristics of nonlinear dynamic friction. According to the velocity curve, the overall velocity tracing accuracy of the robot is also higher than that at low speed. Therefore, the proposed dynamic compensation method can be extended to improve the stability and tracking performance for other type of robots.

Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.