# Designation and Multiobjective Optimization of a New Six-DOF Haptic Device Based on Genetic Algorithm

Authors: Baoyu Shi; Hongtao Wu
Journal: Wireless Communications and Mobile Computing
 (2021)

---

## Abstract

This paper focuses on designation and multiobjective optimization of a new haptic device. A new six-DOF hybrid haptic device has been designed combining series mechanism 3PRPaR and parallel mechanism ROBO_003.3PRPaR which is the positioning mechanism and ROBO_003 which is the directional mechanism. Screw theory was carried out to analyze the DOF of ROBO_003, aiming at the requirement performances of working space, stiffness, reverse drive, and output capacity of haptic device. Forward kinematics and inverse kinematics of the positioning mechanism were analyzed. To ensure the force transmission performance and compact structure, the evaluating indicator is dexterity and volume. Generalized coordinates as constraint condition, length of moving rod, and dexterity were chosen as the objective function. The parameters were optimized based on multiobjective optimization of the genetic algorithm. The paretooptimal solution set of the objective function was received. The results showed that the most suitable parameters are 0.502 rad. The optimized dexterity is 69% higher than the unoptimized, and the rod length is reduced by 6%.The research of this paper provides a theoretical basis for the designation of new haptic device.


---

## Body

## 1. Introduction

With the development of virtual reality technology, the importance of haptic device is becoming increasingly prominent. As an indispensable part of the interactive exploration system, haptic device can realize the movement and force transmission information between human and virtual environment, reproduce the contact between human and environment, and enhance the telepresence of virtual reality. Haptic devices are divided into finger type, arm type, and whole body type. Immersion has developed a finger type haptic device called CyberGrasp [1], which can make operators feel the size and shape of virtual object. Its disadvantage is that it can only control the movement of operators’ whole finger, but cannot give force feedback to a single joint of the finger. Rutgers University of America has developed a built-in multifinger force feedback glove Master II-ND [2], which can measure the posture and feedback force of the haptic device simultaneously. It has the advantages of small mass, large dynamic range, and small friction. The disadvantage of Master II-ND is that it is installed on the palm, the palm cannot hold tightly, and the finger’s range of motion is limited. Among the arm type force feedback devices, the phantom force feedback device of sensable company is more representative [3]. The whole body force feedback equipment is usually bulky and inconvenient to operate. Haptic devices are divided into series structure and parallel structure according to robot configuration [4]. The series structure force feedback equipment which is widely used is phantom series products of sensable technology company [5]. Haption Virtuose introduced by literature [6] was developed by the Haption company of France. It is also a 6-DOF series force feedback device, which can be operated in the range of human multi arm movement. Parallel mechanism has a closed chain structure, which has higher structural stiffness, greater bearing capacity, and higher positioning accuracy [7]. The characteristics make parallel mechanisms meet the requirements of haptic device. The 6-DOF hybrid haptic device introduced in literature [8] was developed by Niigata University of Japan, and 3-DOF in moving directions is provided by delta mechanism. It has a compact structure and is suitable for operation on desktop platform [8], but it has a small stiffness. If the feedback force is large, it may cause small deformation of the mechanism which can affect the accuracy of position detection. Therefore, it is necessary to adopt the mechanism with higher structural stiffness. At present, most of the research on haptic devices are based on surgical or medical robot, which has a certain specificity in the designation of mechanism (including the number of degrees of freedom, the shape, and size of workspace). For example, high feedback force or torque is not needed in surgical robot. Therefore, the output capacity may be insufficient in other applications. Therefore, the existing structure and hardware need to be improved to enhance its versatility.Better performance can be obtained by optimizing the parameters of mechanism. The construction of reasonable kinematic index is an important basis to evaluate the kinematic performance of mechanism. In the structural optimization design of parallel mechanism, Jacobian matrix condition number is used as the evaluation index of mechanism design performance in literature [9]. The singular value of Jacobian matrix is used as the performance index to optimize the structure of parallel mechanism in literature [10]. Carbon used the binary method to optimize the workspace volume of reachable mechanism in literature [11].The dimension synthesis of mechanism is based on the performance index to optimize the size and other parameters of the mechanism [12]. Many literatures focus on the intelligent optimization algorithm. Zhang et al. used the differential evolution algorithm to optimize the stiffness and reachable workspace of 3-UPU in literature [13]. Cui and others optimized the workspace and dexterity of 3-ups-s based on Isight in literature [14]. ZHANG et al. used the genetic algorithm to optimize the global condition number and workspace of spherical mechanism [15].This paper is aimed at the problems of the existing haptic devices, such as low structural stiffness, limited output capacity, and application simplification, and a 6-DOF haptic device with large output capability was designed for general operation object. The 6-DOF haptic device satisfies high structural stiffness and is a universal haptic device. Combined with screw theory, Based on the analysis of its positioning mechanism and orientation mechanism, taking the generalized coordinatesθand β as constraints, the length ρ and dexterity QD of the moving rod as objective functions, the size, volume, and dexterity of the mechanism were optimized in the multiobjective optimization method based on the genetic algorithm. The minimum workspace and the optimal dexterity were obtained, which provide a theoretical basis for the design of universal haptic device with better performance.The following content is structured as follows. The designation and kinematic analysis of the haptic device is in Section2. Section 3 discusses the workspace analysis of the haptic device. Section 4 is structural parameter optimization of the haptic device. Section 5 is multiobjective optimization of GA. Section 6 is results and discussion. Section 7 concludes this paper.


## 2. Designation and Kinematic Analysis of the Haptic Device

### 2.1. Mechanism Designation

The structure of haptic device can be divided into series connection, parallel connection, and hybrid connection. Series mechanism belongs to the open chain structure, and the parallel mechanism contains at least two independent kinematic chains between the end and the base. Series mechanism belongs to the open chain structure, and there are at least two independent kinematic chains between the end and the base of the parallel mechanism. The hybrid mechanism includes series mechanism and parallel mechanism. The advantage of series mechanism is that it has a large working space, but its actuator is installed at the joint of each connecting rod. It increases the motion inertia of the mechanism, which is not conducive to the reverse drive and affects the output ability of the series mechanism. The actuator of the parallel mechanism is installed on the base, which reduces the motion inertia. It has a high structural stiffness and can transfer large force or torque. However, the parallel mechanism is not suitable for force feedback because of its small workspace and insufficient speed output capability. In view of the advantages and disadvantages of the series mechanism and parallel mechanism, a new type of haptic device was designed by using the hybrid structure of the series and parallel mechanism, which is convenient for real-time control of haptic device.


#### 2.1.1. The Designation of the Positioning Mechanism

Generally speaking, the workspace of haptic device is smaller than that of the virtual simulation environment. When mapping the operation space, if the workspace is inconsistent, there will be problems in the mapping of displacement scale factor and reachable workspace. To eliminate these situations, the workspace of the haptic device should be a regular cube space, which can correspond to the real or virtual three-dimensional space. Meanwhile, the structure of the whole mechanism should be kept as compact as possible. In order to achieve fast motion, the moving parts of the haptic device need a small inertia.If there is no force sensor at the end effector, the haptic device needs to be structurally isotropic [16], which ensure good control of forces in all directions. On the contrary, if a force sensor is used, the inertia at the end effector and the size of the haptic device will be increased, and the cost will also be increased. Therefore, it is better to choose the mechanism with isotropic configuration in the workspace. A fully isotropic 3-DOF mechanism was described in literature [17]. The movement direction of the branch chain is perpendicular to the length direction of the rod, which needs to increase the structure size to increase the bending capacity and the motion inertia; so, it is not suitable for haptic device. An improved translational mechanism 3PRPaR was described in literature [18] as shown in Figure 1.Figure 1
3PRPaR.The direction of the moving pair of the mechanism is consistent with the direction of the rod length, and the force transmission direction is also along the direction of the rod length; so, it has high stiffness [19], which can achieve large output of feedback force. At the same time, because the mobile joints are installed perpendicular to each other, they have an isotropic configuration in their workspace. In addition, the workspace of the mechanism is approximately a regular cube; therefore, the mechanism is suitable for the positioning mechanism of haptic device [20]. The branch chain structure of 3PRPaR and schematic diagram of three branch chain combination structures are shown as Figures 2 and 3.Figure 2
PRPaR.Figure 3
3PRPaR.Three branch chains of 3PRPaR can move alongX, Y, and Z directions; so, the positioning mechanism has three degrees of freedom.


#### 2.1.2. The Designation of the Directional Mechanism and the Analysis of DOF

(1)
A 6-DOF parallel mechanism Triplanar (Figure4) was described in literature [21]. Triplanar is driven by three planar motors with two degrees of freedom. The design concept of ROBO_003 comes from Triplanar but has its own obvious characteristics (Figure 5). Three planar motors with two degrees of freedom in Triplanar are replaced by three turntables rotating around the fixed axis in ROBO_003. It is equivalent to limiting the two degrees of freedom motion of three planar motors to a single degree of freedom motion along a fixed circle(2)
Analysis of freedom degree of directional mechanism based on screw theoryFigure 4
Triplanar.Figure 5
ROBO_003.Plu¨ckerHomogeneous coordinates [22]:
(1)S=ss0=sr×s+hs,S=lmnpqrT,l2+m2+n2=1,lp+mq+nr=0,whenS is a motion screw: S=ωv, to get a branch chain P1G1Q1 in ROBO_003. Motion screw of each joint in the branch chain is shown in Figure 6.Figure 6
Screw of the branch chain in ROBO_003.In the branch chainP1G1Q1, P1 is the rotating pair, G1 is the spherical pair, and Q1 is the rotating pair. Five screws in the branch chain are shown in Figure 6 as follows:
(2)s1=001;000,s2=100;00‐b,s3=010;00a,s4=001;b‐a0,s5=def;ghi.The motion screw ofP1G1Q1 in the parallel robot ROBO_003 can be set as S=s1s2s3s4s5. It can be expressed as the matrix:
(3)S=0100d0010e1001f000bg000−ah0−ba0i,STESr=022.ST is the transpose matrix of S, E is the unit matrix, and Sr is the unit force screw of terminal constraint.LetR=STE, and the problem of solving Sr can be transformed into the problem of solving the following homogeneous linear equation.
(4)Rx=0;x=Sr,Sr=x1x2x3;x4x5x6.It was calculated by MATLAB:(5)x1=1;x2=b/a;x3=g+b/ahae−bd−i;x4=bg+b2/ahae−bd−i;x5=−ag+bhae−bd−i;x6=0.According to equation (5), it can be concluded that
(6)x1⋅x4+x2⋅x5+x3⋅x6=0.TheSr is a line vector. The constraints of the other two branched chains of ROBO_003 on the moving platform are force vectors. It is proved that the three force vectors are linear independent. Therefore, the three branches restrict the three degrees of freedom of the moving platform. For the spatial mechanism, only three degrees of freedom of rotation are unconstrained. Therefore, ROBO_003 is a pure rotating mechanism with three degrees of freedom. Its moving platform can rotate along the x, y, and z axis, respectively.The movement and rotation of the haptic device are decoupled. The assembly diagram of the haptic device is shown in Figure7.Figure 7
Assembly diagram of the haptic device.


### 2.2. Kinematic Analysis of the Haptic Device

It is necessary to analyze the kinematics of the haptic device and establish the mapping relationship between the joint space and the operation space of the haptic device. Deduce the position of the end effector in the space by detecting the position change of its joint and then transfer the position information of the haptic device to the slave and carry out the inverse kinematics. Control the movement of the slave joints to make the slave reach the target position.


#### 2.2.1. Kinematic Analysis of the Positioning Mechanism

(1) Inverse Kinematic Analysis. 3PRPaR consists of three branch chains (Figure 3), and the prismatic pair of the three branch chains is perpendicular to each other. Point O is defined as the origin of fixed coordinate system and P as the origin of the moving platform coordinate system. Point P coincides with point O at initial moment. The X axis, Y axis, and Z axis point to the origin along the prismatic pair, respectively. As shown in Figure 7, the three branch chains are spatially symmetric. So, we only need to analyze the kinematics of a branch chain. In order to simplify the calculation, the parallelogram joint in each branch chain is replaced by an equivalent branch chain composed of a connecting rod and two universal joints. Ai1 (i=1,2,3) is driving slider of linear motion pair as shown in Figure 8.Figure 8
Equivalent kinematic chain of the single chain.Given generalized coordinatesβi1, βi2, θi1, θi2, and βi1 are the angles that rotate around the axis zi1, βi2 is the angle that rotate around the axis zi5, θi1 is the angle that rotate around the axis zi3, and θi2 is the angle that rotate around the axis zi4.(i is the ith branch chain): βi1=-βi2, θi1=-θi2. The end position of the first branch chain can be expressed as
(7)x1y1z1=q11+d4cosβ1cosθ1+d6−d4sinβ1−d4cosβ1sinθ1.The end positions of branched chain 2 and branched chain 3 can also be expressed as(8)x2y2z2=d4cosβ2sinθ2q21+d4cosβ2cosθ2+d6d4sinβ2.(9)x3y3z3=d4sinβ3d4cosβ3sinθ3q31+d4cosβ3cosθ3+d6.The ends of the three branch chains are connected together, so that(10)xyz=x1y1z1=x2y2z2=x3y3z3.Equation (10) square is as follows:
(11)x2y2z2=q11+d4cosβ1cosθ1+d62−d4sinβ12−d4cosβ1sinθ12.Each line of equation (11) is added and simplified:
(12)x−q11−d62+y2+z2−d42=0.Similarly, it was added and simplified the lines of equation (8) and (9), respectively:
(13)x2+y−q21−d62+z2−d42=0;x2+y2+z−q31−d62−d42=0.To solveqi1, expand each lines of equation (13) and simplified:
(14)qi12+ci1qi1+ci0=0,i=1,2,3;c11=2d6−2x;c10=−d42+d62−2xd6+x2+y2+z2;c21=2d6−2yc20=−d42+d62−2yd6+x2+y2+z2c31=2d6−2zc30=−d42+d62−2zd6+x2+y2+z2.Solveqi12+ci1qi1+ci0=0,i=1,2,3 and get
(15)qi1=−ci1±ci12−4ci02i=1,2,3.“±” means two different working patterns, as shown in Figures 9 and 10.Figure 9
Working pattern 1.Figure 10
Working pattern 2.In Figure9, branch chain 1 has crossed the singular boundary point, and the travel has exceeded the limit displacement range of the guide rail. Therefore, the displacement can be used to determine the unique solution of the inverse kinematic solution. After the variables of each driving joint are calculated, the variable value of the passive joint can be obtained according to the geometric relationship.
(16)θ1=atan2−z,x−q11−d6;θ2=atan2x,y−q21−d6;θ3=atan2y,z−q31−d6,β1=atan2−y,x−q11−d6cosθ1;β2=atan2z,y−q21−d6cosθ2β3=atan2x,z−q31−d6cosθ3.(2) Forward Kinematic Analysis. Equation (10) can be written as follows:
(17)x2+y2+z2+2c1x+2c2z+c12+c22−d42x2+y2+z2+2c3x+2c4y+2c5z+c32+c42+c52−d42x2+y2+z2+2c6x+2c7y+2c8z+c62+c72+c82−d42=000,(18)c1=−d6,c2=−q11,c3=−q21,c4=−d6,c4=−d6,c5=0,c6=−d6,c7=−q31,c8=0.Solve equations (17), and we get the forward kinematic model of 3PRRaR:
(19)x=c9z+c10,y=c11z+c12,z=−c14±c142−4c13c152c13,c9=−c2−c5c7c1c4+c1c7−c3c7−c4c6;c10=−c12+c22−c32−c42−c52c4c1c4+c1c7−c3c7−c4c6;c11=c2−c8c1−c6c1c4+c1c7−c3c7−c4c6,c12=−c12+c22−c62−c72−c82c1−c3c1c4+c1c7−c3c7−c4c6;c13=c92+c112+1;c14=2c9c10+c11c12+c1c9−c2,c15=c102+c122+2c1c10+c12+c22−d42.


## 2.1. Mechanism Designation

The structure of haptic device can be divided into series connection, parallel connection, and hybrid connection. Series mechanism belongs to the open chain structure, and the parallel mechanism contains at least two independent kinematic chains between the end and the base. Series mechanism belongs to the open chain structure, and there are at least two independent kinematic chains between the end and the base of the parallel mechanism. The hybrid mechanism includes series mechanism and parallel mechanism. The advantage of series mechanism is that it has a large working space, but its actuator is installed at the joint of each connecting rod. It increases the motion inertia of the mechanism, which is not conducive to the reverse drive and affects the output ability of the series mechanism. The actuator of the parallel mechanism is installed on the base, which reduces the motion inertia. It has a high structural stiffness and can transfer large force or torque. However, the parallel mechanism is not suitable for force feedback because of its small workspace and insufficient speed output capability. In view of the advantages and disadvantages of the series mechanism and parallel mechanism, a new type of haptic device was designed by using the hybrid structure of the series and parallel mechanism, which is convenient for real-time control of haptic device.


### 2.1.1. The Designation of the Positioning Mechanism

Generally speaking, the workspace of haptic device is smaller than that of the virtual simulation environment. When mapping the operation space, if the workspace is inconsistent, there will be problems in the mapping of displacement scale factor and reachable workspace. To eliminate these situations, the workspace of the haptic device should be a regular cube space, which can correspond to the real or virtual three-dimensional space. Meanwhile, the structure of the whole mechanism should be kept as compact as possible. In order to achieve fast motion, the moving parts of the haptic device need a small inertia.If there is no force sensor at the end effector, the haptic device needs to be structurally isotropic [16], which ensure good control of forces in all directions. On the contrary, if a force sensor is used, the inertia at the end effector and the size of the haptic device will be increased, and the cost will also be increased. Therefore, it is better to choose the mechanism with isotropic configuration in the workspace. A fully isotropic 3-DOF mechanism was described in literature [17]. The movement direction of the branch chain is perpendicular to the length direction of the rod, which needs to increase the structure size to increase the bending capacity and the motion inertia; so, it is not suitable for haptic device. An improved translational mechanism 3PRPaR was described in literature [18] as shown in Figure 1.Figure 1
3PRPaR.The direction of the moving pair of the mechanism is consistent with the direction of the rod length, and the force transmission direction is also along the direction of the rod length; so, it has high stiffness [19], which can achieve large output of feedback force. At the same time, because the mobile joints are installed perpendicular to each other, they have an isotropic configuration in their workspace. In addition, the workspace of the mechanism is approximately a regular cube; therefore, the mechanism is suitable for the positioning mechanism of haptic device [20]. The branch chain structure of 3PRPaR and schematic diagram of three branch chain combination structures are shown as Figures 2 and 3.Figure 2
PRPaR.Figure 3
3PRPaR.Three branch chains of 3PRPaR can move alongX, Y, and Z directions; so, the positioning mechanism has three degrees of freedom.


### 2.1.2. The Designation of the Directional Mechanism and the Analysis of DOF

(1)
A 6-DOF parallel mechanism Triplanar (Figure4) was described in literature [21]. Triplanar is driven by three planar motors with two degrees of freedom. The design concept of ROBO_003 comes from Triplanar but has its own obvious characteristics (Figure 5). Three planar motors with two degrees of freedom in Triplanar are replaced by three turntables rotating around the fixed axis in ROBO_003. It is equivalent to limiting the two degrees of freedom motion of three planar motors to a single degree of freedom motion along a fixed circle(2)
Analysis of freedom degree of directional mechanism based on screw theoryFigure 4
Triplanar.Figure 5
ROBO_003.Plu¨ckerHomogeneous coordinates [22]:
(1)S=ss0=sr×s+hs,S=lmnpqrT,l2+m2+n2=1,lp+mq+nr=0,whenS is a motion screw: S=ωv, to get a branch chain P1G1Q1 in ROBO_003. Motion screw of each joint in the branch chain is shown in Figure 6.Figure 6
Screw of the branch chain in ROBO_003.In the branch chainP1G1Q1, P1 is the rotating pair, G1 is the spherical pair, and Q1 is the rotating pair. Five screws in the branch chain are shown in Figure 6 as follows:
(2)s1=001;000,s2=100;00‐b,s3=010;00a,s4=001;b‐a0,s5=def;ghi.The motion screw ofP1G1Q1 in the parallel robot ROBO_003 can be set as S=s1s2s3s4s5. It can be expressed as the matrix:
(3)S=0100d0010e1001f000bg000−ah0−ba0i,STESr=022.ST is the transpose matrix of S, E is the unit matrix, and Sr is the unit force screw of terminal constraint.LetR=STE, and the problem of solving Sr can be transformed into the problem of solving the following homogeneous linear equation.
(4)Rx=0;x=Sr,Sr=x1x2x3;x4x5x6.It was calculated by MATLAB:(5)x1=1;x2=b/a;x3=g+b/ahae−bd−i;x4=bg+b2/ahae−bd−i;x5=−ag+bhae−bd−i;x6=0.According to equation (5), it can be concluded that
(6)x1⋅x4+x2⋅x5+x3⋅x6=0.TheSr is a line vector. The constraints of the other two branched chains of ROBO_003 on the moving platform are force vectors. It is proved that the three force vectors are linear independent. Therefore, the three branches restrict the three degrees of freedom of the moving platform. For the spatial mechanism, only three degrees of freedom of rotation are unconstrained. Therefore, ROBO_003 is a pure rotating mechanism with three degrees of freedom. Its moving platform can rotate along the x, y, and z axis, respectively.The movement and rotation of the haptic device are decoupled. The assembly diagram of the haptic device is shown in Figure7.Figure 7
Assembly diagram of the haptic device.


## 2.1.1. The Designation of the Positioning Mechanism

Generally speaking, the workspace of haptic device is smaller than that of the virtual simulation environment. When mapping the operation space, if the workspace is inconsistent, there will be problems in the mapping of displacement scale factor and reachable workspace. To eliminate these situations, the workspace of the haptic device should be a regular cube space, which can correspond to the real or virtual three-dimensional space. Meanwhile, the structure of the whole mechanism should be kept as compact as possible. In order to achieve fast motion, the moving parts of the haptic device need a small inertia.If there is no force sensor at the end effector, the haptic device needs to be structurally isotropic [16], which ensure good control of forces in all directions. On the contrary, if a force sensor is used, the inertia at the end effector and the size of the haptic device will be increased, and the cost will also be increased. Therefore, it is better to choose the mechanism with isotropic configuration in the workspace. A fully isotropic 3-DOF mechanism was described in literature [17]. The movement direction of the branch chain is perpendicular to the length direction of the rod, which needs to increase the structure size to increase the bending capacity and the motion inertia; so, it is not suitable for haptic device. An improved translational mechanism 3PRPaR was described in literature [18] as shown in Figure 1.Figure 1
3PRPaR.The direction of the moving pair of the mechanism is consistent with the direction of the rod length, and the force transmission direction is also along the direction of the rod length; so, it has high stiffness [19], which can achieve large output of feedback force. At the same time, because the mobile joints are installed perpendicular to each other, they have an isotropic configuration in their workspace. In addition, the workspace of the mechanism is approximately a regular cube; therefore, the mechanism is suitable for the positioning mechanism of haptic device [20]. The branch chain structure of 3PRPaR and schematic diagram of three branch chain combination structures are shown as Figures 2 and 3.Figure 2
PRPaR.Figure 3
3PRPaR.Three branch chains of 3PRPaR can move alongX, Y, and Z directions; so, the positioning mechanism has three degrees of freedom.


## 2.1.2. The Designation of the Directional Mechanism and the Analysis of DOF

(1)
A 6-DOF parallel mechanism Triplanar (Figure4) was described in literature [21]. Triplanar is driven by three planar motors with two degrees of freedom. The design concept of ROBO_003 comes from Triplanar but has its own obvious characteristics (Figure 5). Three planar motors with two degrees of freedom in Triplanar are replaced by three turntables rotating around the fixed axis in ROBO_003. It is equivalent to limiting the two degrees of freedom motion of three planar motors to a single degree of freedom motion along a fixed circle(2)
Analysis of freedom degree of directional mechanism based on screw theoryFigure 4
Triplanar.Figure 5
ROBO_003.Plu¨ckerHomogeneous coordinates [22]:
(1)S=ss0=sr×s+hs,S=lmnpqrT,l2+m2+n2=1,lp+mq+nr=0,whenS is a motion screw: S=ωv, to get a branch chain P1G1Q1 in ROBO_003. Motion screw of each joint in the branch chain is shown in Figure 6.Figure 6
Screw of the branch chain in ROBO_003.In the branch chainP1G1Q1, P1 is the rotating pair, G1 is the spherical pair, and Q1 is the rotating pair. Five screws in the branch chain are shown in Figure 6 as follows:
(2)s1=001;000,s2=100;00‐b,s3=010;00a,s4=001;b‐a0,s5=def;ghi.The motion screw ofP1G1Q1 in the parallel robot ROBO_003 can be set as S=s1s2s3s4s5. It can be expressed as the matrix:
(3)S=0100d0010e1001f000bg000−ah0−ba0i,STESr=022.ST is the transpose matrix of S, E is the unit matrix, and Sr is the unit force screw of terminal constraint.LetR=STE, and the problem of solving Sr can be transformed into the problem of solving the following homogeneous linear equation.
(4)Rx=0;x=Sr,Sr=x1x2x3;x4x5x6.It was calculated by MATLAB:(5)x1=1;x2=b/a;x3=g+b/ahae−bd−i;x4=bg+b2/ahae−bd−i;x5=−ag+bhae−bd−i;x6=0.According to equation (5), it can be concluded that
(6)x1⋅x4+x2⋅x5+x3⋅x6=0.TheSr is a line vector. The constraints of the other two branched chains of ROBO_003 on the moving platform are force vectors. It is proved that the three force vectors are linear independent. Therefore, the three branches restrict the three degrees of freedom of the moving platform. For the spatial mechanism, only three degrees of freedom of rotation are unconstrained. Therefore, ROBO_003 is a pure rotating mechanism with three degrees of freedom. Its moving platform can rotate along the x, y, and z axis, respectively.The movement and rotation of the haptic device are decoupled. The assembly diagram of the haptic device is shown in Figure7.Figure 7
Assembly diagram of the haptic device.


## 2.2. Kinematic Analysis of the Haptic Device

It is necessary to analyze the kinematics of the haptic device and establish the mapping relationship between the joint space and the operation space of the haptic device. Deduce the position of the end effector in the space by detecting the position change of its joint and then transfer the position information of the haptic device to the slave and carry out the inverse kinematics. Control the movement of the slave joints to make the slave reach the target position.


### 2.2.1. Kinematic Analysis of the Positioning Mechanism

(1) Inverse Kinematic Analysis. 3PRPaR consists of three branch chains (Figure 3), and the prismatic pair of the three branch chains is perpendicular to each other. Point O is defined as the origin of fixed coordinate system and P as the origin of the moving platform coordinate system. Point P coincides with point O at initial moment. The X axis, Y axis, and Z axis point to the origin along the prismatic pair, respectively. As shown in Figure 7, the three branch chains are spatially symmetric. So, we only need to analyze the kinematics of a branch chain. In order to simplify the calculation, the parallelogram joint in each branch chain is replaced by an equivalent branch chain composed of a connecting rod and two universal joints. Ai1 (i=1,2,3) is driving slider of linear motion pair as shown in Figure 8.Figure 8
Equivalent kinematic chain of the single chain.Given generalized coordinatesβi1, βi2, θi1, θi2, and βi1 are the angles that rotate around the axis zi1, βi2 is the angle that rotate around the axis zi5, θi1 is the angle that rotate around the axis zi3, and θi2 is the angle that rotate around the axis zi4.(i is the ith branch chain): βi1=-βi2, θi1=-θi2. The end position of the first branch chain can be expressed as
(7)x1y1z1=q11+d4cosβ1cosθ1+d6−d4sinβ1−d4cosβ1sinθ1.The end positions of branched chain 2 and branched chain 3 can also be expressed as(8)x2y2z2=d4cosβ2sinθ2q21+d4cosβ2cosθ2+d6d4sinβ2.(9)x3y3z3=d4sinβ3d4cosβ3sinθ3q31+d4cosβ3cosθ3+d6.The ends of the three branch chains are connected together, so that(10)xyz=x1y1z1=x2y2z2=x3y3z3.Equation (10) square is as follows:
(11)x2y2z2=q11+d4cosβ1cosθ1+d62−d4sinβ12−d4cosβ1sinθ12.Each line of equation (11) is added and simplified:
(12)x−q11−d62+y2+z2−d42=0.Similarly, it was added and simplified the lines of equation (8) and (9), respectively:
(13)x2+y−q21−d62+z2−d42=0;x2+y2+z−q31−d62−d42=0.To solveqi1, expand each lines of equation (13) and simplified:
(14)qi12+ci1qi1+ci0=0,i=1,2,3;c11=2d6−2x;c10=−d42+d62−2xd6+x2+y2+z2;c21=2d6−2yc20=−d42+d62−2yd6+x2+y2+z2c31=2d6−2zc30=−d42+d62−2zd6+x2+y2+z2.Solveqi12+ci1qi1+ci0=0,i=1,2,3 and get
(15)qi1=−ci1±ci12−4ci02i=1,2,3.“±” means two different working patterns, as shown in Figures 9 and 10.Figure 9
Working pattern 1.Figure 10
Working pattern 2.In Figure9, branch chain 1 has crossed the singular boundary point, and the travel has exceeded the limit displacement range of the guide rail. Therefore, the displacement can be used to determine the unique solution of the inverse kinematic solution. After the variables of each driving joint are calculated, the variable value of the passive joint can be obtained according to the geometric relationship.
(16)θ1=atan2−z,x−q11−d6;θ2=atan2x,y−q21−d6;θ3=atan2y,z−q31−d6,β1=atan2−y,x−q11−d6cosθ1;β2=atan2z,y−q21−d6cosθ2β3=atan2x,z−q31−d6cosθ3.(2) Forward Kinematic Analysis. Equation (10) can be written as follows:
(17)x2+y2+z2+2c1x+2c2z+c12+c22−d42x2+y2+z2+2c3x+2c4y+2c5z+c32+c42+c52−d42x2+y2+z2+2c6x+2c7y+2c8z+c62+c72+c82−d42=000,(18)c1=−d6,c2=−q11,c3=−q21,c4=−d6,c4=−d6,c5=0,c6=−d6,c7=−q31,c8=0.Solve equations (17), and we get the forward kinematic model of 3PRRaR:
(19)x=c9z+c10,y=c11z+c12,z=−c14±c142−4c13c152c13,c9=−c2−c5c7c1c4+c1c7−c3c7−c4c6;c10=−c12+c22−c32−c42−c52c4c1c4+c1c7−c3c7−c4c6;c11=c2−c8c1−c6c1c4+c1c7−c3c7−c4c6,c12=−c12+c22−c62−c72−c82c1−c3c1c4+c1c7−c3c7−c4c6;c13=c92+c112+1;c14=2c9c10+c11c12+c1c9−c2,c15=c102+c122+2c1c10+c12+c22−d42.


## 2.2.1. Kinematic Analysis of the Positioning Mechanism

(1) Inverse Kinematic Analysis. 3PRPaR consists of three branch chains (Figure 3), and the prismatic pair of the three branch chains is perpendicular to each other. Point O is defined as the origin of fixed coordinate system and P as the origin of the moving platform coordinate system. Point P coincides with point O at initial moment. The X axis, Y axis, and Z axis point to the origin along the prismatic pair, respectively. As shown in Figure 7, the three branch chains are spatially symmetric. So, we only need to analyze the kinematics of a branch chain. In order to simplify the calculation, the parallelogram joint in each branch chain is replaced by an equivalent branch chain composed of a connecting rod and two universal joints. Ai1 (i=1,2,3) is driving slider of linear motion pair as shown in Figure 8.Figure 8
Equivalent kinematic chain of the single chain.Given generalized coordinatesβi1, βi2, θi1, θi2, and βi1 are the angles that rotate around the axis zi1, βi2 is the angle that rotate around the axis zi5, θi1 is the angle that rotate around the axis zi3, and θi2 is the angle that rotate around the axis zi4.(i is the ith branch chain): βi1=-βi2, θi1=-θi2. The end position of the first branch chain can be expressed as
(7)x1y1z1=q11+d4cosβ1cosθ1+d6−d4sinβ1−d4cosβ1sinθ1.The end positions of branched chain 2 and branched chain 3 can also be expressed as(8)x2y2z2=d4cosβ2sinθ2q21+d4cosβ2cosθ2+d6d4sinβ2.(9)x3y3z3=d4sinβ3d4cosβ3sinθ3q31+d4cosβ3cosθ3+d6.The ends of the three branch chains are connected together, so that(10)xyz=x1y1z1=x2y2z2=x3y3z3.Equation (10) square is as follows:
(11)x2y2z2=q11+d4cosβ1cosθ1+d62−d4sinβ12−d4cosβ1sinθ12.Each line of equation (11) is added and simplified:
(12)x−q11−d62+y2+z2−d42=0.Similarly, it was added and simplified the lines of equation (8) and (9), respectively:
(13)x2+y−q21−d62+z2−d42=0;x2+y2+z−q31−d62−d42=0.To solveqi1, expand each lines of equation (13) and simplified:
(14)qi12+ci1qi1+ci0=0,i=1,2,3;c11=2d6−2x;c10=−d42+d62−2xd6+x2+y2+z2;c21=2d6−2yc20=−d42+d62−2yd6+x2+y2+z2c31=2d6−2zc30=−d42+d62−2zd6+x2+y2+z2.Solveqi12+ci1qi1+ci0=0,i=1,2,3 and get
(15)qi1=−ci1±ci12−4ci02i=1,2,3.“±” means two different working patterns, as shown in Figures 9 and 10.Figure 9
Working pattern 1.Figure 10
Working pattern 2.In Figure9, branch chain 1 has crossed the singular boundary point, and the travel has exceeded the limit displacement range of the guide rail. Therefore, the displacement can be used to determine the unique solution of the inverse kinematic solution. After the variables of each driving joint are calculated, the variable value of the passive joint can be obtained according to the geometric relationship.
(16)θ1=atan2−z,x−q11−d6;θ2=atan2x,y−q21−d6;θ3=atan2y,z−q31−d6,β1=atan2−y,x−q11−d6cosθ1;β2=atan2z,y−q21−d6cosθ2β3=atan2x,z−q31−d6cosθ3.(2) Forward Kinematic Analysis. Equation (10) can be written as follows:
(17)x2+y2+z2+2c1x+2c2z+c12+c22−d42x2+y2+z2+2c3x+2c4y+2c5z+c32+c42+c52−d42x2+y2+z2+2c6x+2c7y+2c8z+c62+c72+c82−d42=000,(18)c1=−d6,c2=−q11,c3=−q21,c4=−d6,c4=−d6,c5=0,c6=−d6,c7=−q31,c8=0.Solve equations (17), and we get the forward kinematic model of 3PRRaR:
(19)x=c9z+c10,y=c11z+c12,z=−c14±c142−4c13c152c13,c9=−c2−c5c7c1c4+c1c7−c3c7−c4c6;c10=−c12+c22−c32−c42−c52c4c1c4+c1c7−c3c7−c4c6;c11=c2−c8c1−c6c1c4+c1c7−c3c7−c4c6,c12=−c12+c22−c62−c72−c82c1−c3c1c4+c1c7−c3c7−c4c6;c13=c92+c112+1;c14=2c9c10+c11c12+c1c9−c2,c15=c102+c122+2c1c10+c12+c22−d42.


## 3. Workspace Analysis of the Haptic Device

The workspace of haptic device should be cube, which is convenient to establish the spatial mapping relationship between master hand and slave hand in the teleoperation system.


### 3.1. Monte Carlo Method

The Monte Carlo method is a numerical calculation method based on statistical theory and probability theory [23]. The basic idea is to use random number to solve the problem with definite law. When it is used to solve the problem of robot workspace, the joint variables can be used as input, and the space position of the end of the manipulator can be obtained through the forward kinematic calculation of the robot. When the number of joint variables is enough, the set of corresponding end space positions will approximately constitute the workspace of the robot.The steps are as follows:(1)
Determine the range of motion of each joint(2)
The forward kinematics of the mobile mechanism is analyzed(3)
By randomly selecting an input for each joint, the random variables of each joint can be obtained by the following formula:θi=θi−min+θi−max−θi−min×randθi−max is the maximum input angle of the ith joint.θi−min is the minimum input angle of the ith joint.Rand is the function that generates random numbers(4)
The randomly generated joint variables are brought into the kinematic equation to calculate the position of the operation end in the operation space(5)
Repeat the above steps to get the set of spatial positions of the operation end, which is the approximate solution of the workspace. The more groups of joint variables are taken, the closer the result is to the actual workspaceFollow the above steps. The workspace of the positioning mechanism is calculated by MATLAB as shown in Figures11–14.Figure 11
Workspace of the haptic device.Figure 12
X-Y view.Figure 13
X=0 slice diagram.Figure 14
Y=0 slice diagram.It contains a square in Figure14.So, in the workspace of 3PRPaR, it contains a cube. It meets the design requirements of the haptic device.


## 3.1. Monte Carlo Method

The Monte Carlo method is a numerical calculation method based on statistical theory and probability theory [23]. The basic idea is to use random number to solve the problem with definite law. When it is used to solve the problem of robot workspace, the joint variables can be used as input, and the space position of the end of the manipulator can be obtained through the forward kinematic calculation of the robot. When the number of joint variables is enough, the set of corresponding end space positions will approximately constitute the workspace of the robot.The steps are as follows:(1)
Determine the range of motion of each joint(2)
The forward kinematics of the mobile mechanism is analyzed(3)
By randomly selecting an input for each joint, the random variables of each joint can be obtained by the following formula:θi=θi−min+θi−max−θi−min×randθi−max is the maximum input angle of the ith joint.θi−min is the minimum input angle of the ith joint.Rand is the function that generates random numbers(4)
The randomly generated joint variables are brought into the kinematic equation to calculate the position of the operation end in the operation space(5)
Repeat the above steps to get the set of spatial positions of the operation end, which is the approximate solution of the workspace. The more groups of joint variables are taken, the closer the result is to the actual workspaceFollow the above steps. The workspace of the positioning mechanism is calculated by MATLAB as shown in Figures11–14.Figure 11
Workspace of the haptic device.Figure 12
X-Y view.Figure 13
X=0 slice diagram.Figure 14
Y=0 slice diagram.It contains a square in Figure14.So, in the workspace of 3PRPaR, it contains a cube. It meets the design requirements of the haptic device.


## 4. Parameter Optimization and Evaluation Index of the Haptic Device

The objective of haptic device optimization is to optimize the size of the mechanism [24]. First of all, the workspace should meet the preset requirements, so that the size of the whole mechanism is as small as possible. At the same time, the transfer coefficient of the mechanism should be limited in a certain range to ensure the stiffness and reverse driving performance of the mechanism. Therefore, it is necessary to select appropriate evaluation indexes. The main evaluation indexes of the haptic device are dexterity and operation space volume.


### 4.1. Dexterity

Dexterity is an important index to measure the transmission accuracy of input and output motion or force of robot [25]. When the robot structure is in or near the singular configuration, the transfer relationship between input motion and output motion will be distorted. The minimum condition number, minimum singular value, and operability based on Jacobian matrix can reflect the degree of motion distortion to a certain extent.The minimum condition number is as follows:(20)kJ=JJ−1.· is the arbitrary matrix norm; when using Frobenius norm, the above formula can be expressed as
(21)kFJ=1ntrJJT1ntrJJT−1.WhenkFJ=1, the mechanism is in isotropic configuration and has the best performance of motion and force transmission. The main axis direction of operating ellipsoid is defined by eigenvalues of JJT−1. Square root of eigenvalue ξ1,ξ2,ξ3 is the length of the spindle. The length of the spindle of operating ellipsoid is directly related to the operating performance of the mechanism. The shorter the length is, the greater the structural stiffness along this direction is. That is, the actuator only needs to output a small force to offset the effect of the external force. Define the velocity transfer coefficients along the three principal axes that are ψ1=1/ξ1,ψ2=1/ξ2,ψ3=1/ξ3 respectively. In order to ensure good force transmission performance and speed transmission performance, ψi needs to be limited to a certain range. Let ψmin≤ψi≤ψmax and ψmin=1/ψmax and get 1/ψmax≤ψi≤ψmax.As shown in Figure15, the workspace of the haptic device contains a cube. Point Q1 and point Q2 are two vertices on the diagonal of the cube. When the haptic device is in singular configuration, point Q1 and pointQ2 are two points closest to the end of the haptic device. According to interval analysis theory [26], if the condition can be satisfied with 1/ψmax≤ψi≤ψmax at points Q1 and Q2, ψi can satisfy this condition in the whole workspace.Figure 15
Singular configuration of the haptic device.As shown in Figures16 and 17, ρii=1,2,3 is the length of the moving rod. When the moving platform of the haptic device is at point Q1, ρ1=ρ2=ρ3=ρmin, when the moving platform of the haptic device is at point Q2, ρ1=ρ2=ρ3=ρmax.Figure 16
Configuration ofQ1.Figure 17
Configuration ofQ2.As shown in Figure18, the length of parallelogram connecting rod is L. The distance from the end of the connecting rod to the center of the moving platform is e. Let ρmin=0, and the coordinate of point P along Z axis can be expressed as Z=−sinβ1L orZ=−sinθ1cosβ2L. In the axis Q1Q2, β1=β2=β3=β, θ1=θ2=θ3=θ. The relationship between β and θ can be expressed as β=−arctansinθ. By substituting the haptic device kinematic equation, we get the following results:
(22)J−1=1−tanθ−tanθ−tanθ11−tanθ−tanθ1,ξ1=2tanθ−1,ξ2=ξ3=tanθ+1ψ1=12tanθ−1,ψ2=ψ3=1tanθ+1.Figure 18
Calculation diagram of the single chain.Solve the equation:(23)1ψmax≤1tanθ+1≤ψmax,1ψmax≤12tanθ−1≤ψmax.QD is the dexterity.
(24)QD=11+tanθ4tanθ−12.


### 4.2. Workspace Volume

As shown in Figure18, a=OA1=OA2=OA3. when point P coincides with point Q1′,
(25)OA2=OQ1′+Q1′C2′+C2′A2.Forρ2=0, C2A2=C2B2=L, Q1′C2=PC2−e, and OQ1′=q1, we get
(26)a=q1+e+L.q1 can be obtained by formula (13); so, a can be obtained bye,L,ψmax. Project A2P=A2B2+B2C2+C2P onto the Y axis gets:
(27)ρ=q1+a+cosθcosβL+e+37,(28)q1=sinβL;a=q1+e+L;L=100sinβ+153.5mm;e=51mm.ρ determines the volume of the workspace. Its length is determined by θ,β.In summary,(29)ρ=q1+a+cosθcosβL+e,(30)QD=11+tanθ4tanθ−12,(31)θ=−arcsintanβ,β∈0π4.


## 4.1. Dexterity

Dexterity is an important index to measure the transmission accuracy of input and output motion or force of robot [25]. When the robot structure is in or near the singular configuration, the transfer relationship between input motion and output motion will be distorted. The minimum condition number, minimum singular value, and operability based on Jacobian matrix can reflect the degree of motion distortion to a certain extent.The minimum condition number is as follows:(20)kJ=JJ−1.· is the arbitrary matrix norm; when using Frobenius norm, the above formula can be expressed as
(21)kFJ=1ntrJJT1ntrJJT−1.WhenkFJ=1, the mechanism is in isotropic configuration and has the best performance of motion and force transmission. The main axis direction of operating ellipsoid is defined by eigenvalues of JJT−1. Square root of eigenvalue ξ1,ξ2,ξ3 is the length of the spindle. The length of the spindle of operating ellipsoid is directly related to the operating performance of the mechanism. The shorter the length is, the greater the structural stiffness along this direction is. That is, the actuator only needs to output a small force to offset the effect of the external force. Define the velocity transfer coefficients along the three principal axes that are ψ1=1/ξ1,ψ2=1/ξ2,ψ3=1/ξ3 respectively. In order to ensure good force transmission performance and speed transmission performance, ψi needs to be limited to a certain range. Let ψmin≤ψi≤ψmax and ψmin=1/ψmax and get 1/ψmax≤ψi≤ψmax.As shown in Figure15, the workspace of the haptic device contains a cube. Point Q1 and point Q2 are two vertices on the diagonal of the cube. When the haptic device is in singular configuration, point Q1 and pointQ2 are two points closest to the end of the haptic device. According to interval analysis theory [26], if the condition can be satisfied with 1/ψmax≤ψi≤ψmax at points Q1 and Q2, ψi can satisfy this condition in the whole workspace.Figure 15
Singular configuration of the haptic device.As shown in Figures16 and 17, ρii=1,2,3 is the length of the moving rod. When the moving platform of the haptic device is at point Q1, ρ1=ρ2=ρ3=ρmin, when the moving platform of the haptic device is at point Q2, ρ1=ρ2=ρ3=ρmax.Figure 16
Configuration ofQ1.Figure 17
Configuration ofQ2.As shown in Figure18, the length of parallelogram connecting rod is L. The distance from the end of the connecting rod to the center of the moving platform is e. Let ρmin=0, and the coordinate of point P along Z axis can be expressed as Z=−sinβ1L orZ=−sinθ1cosβ2L. In the axis Q1Q2, β1=β2=β3=β, θ1=θ2=θ3=θ. The relationship between β and θ can be expressed as β=−arctansinθ. By substituting the haptic device kinematic equation, we get the following results:
(22)J−1=1−tanθ−tanθ−tanθ11−tanθ−tanθ1,ξ1=2tanθ−1,ξ2=ξ3=tanθ+1ψ1=12tanθ−1,ψ2=ψ3=1tanθ+1.Figure 18
Calculation diagram of the single chain.Solve the equation:(23)1ψmax≤1tanθ+1≤ψmax,1ψmax≤12tanθ−1≤ψmax.QD is the dexterity.
(24)QD=11+tanθ4tanθ−12.


## 4.2. Workspace Volume

As shown in Figure18, a=OA1=OA2=OA3. when point P coincides with point Q1′,
(25)OA2=OQ1′+Q1′C2′+C2′A2.Forρ2=0, C2A2=C2B2=L, Q1′C2=PC2−e, and OQ1′=q1, we get
(26)a=q1+e+L.q1 can be obtained by formula (13); so, a can be obtained bye,L,ψmax. Project A2P=A2B2+B2C2+C2P onto the Y axis gets:
(27)ρ=q1+a+cosθcosβL+e+37,(28)q1=sinβL;a=q1+e+L;L=100sinβ+153.5mm;e=51mm.ρ determines the volume of the workspace. Its length is determined by θ,β.In summary,(29)ρ=q1+a+cosθcosβL+e,(30)QD=11+tanθ4tanθ−12,(31)θ=−arcsintanβ,β∈0π4.


## 5. Multiobjective Optimization of the Genetic Algorithm

GA (genetic algorithm) is a kind of evolutionary algorithm. The basic principle is to encode the problem parameters as chromosomes and to exchange the information of chromosomes in the population by means of iterative operation such as selection, crossover, and mutation, so as to generate the chromosomes that meet the optimization objectives [27]. In the genetic algorithm, chromosomes correspond to data or array, which is represented by one-dimensional string structure data. A string of genes is a chromosome which is called individuals. A certain number of individuals form a population, and the number of individuals in population is called population size. Each individual’s adaptability to the environment is called fitness. Flow chart of the genetic algorithm is shown in Figure 19.Figure 19
Flow chart of the genetic algorithm.


### 5.1. Paretooptimal Solution of Multiobjective Optimization

The problem of multiobjective optimization can be described as follows:(32)minf1x,f2x,⋯,fmx,s.t.lb≤x≤ub,Aeq∗x=beq,A∗x≤b.fix is the objective function to be optimized, x (lb≤x≤ub): variables to be optimized.Aeq∗x=beq is the linear equality constraint of x, A∗x≤b: linear inequality constraints.lb≤x≤ub is the range of x.


#### 5.1.1. Gamultiobj Function

The gamultiobj function is an improved multiobjective optimization algorithm based on NSGA (a variant of NSGA). The emergence of the gamultiobj function provides a good way to solve the multiobjective optimization problem in the MATLAB platform. Organization chart of multiobjective optimization function gamultiobj based on GA is shown in Figure20. In this paper, the length ρ and dexterity QD are taken as the objective functions, and the independent variables are θ,β. The best matching value between the volume space suitable for the optimal dexterity is shown in Figures 21–23.Figure 20
Organization chart of the gamultiobj function.Figure 21
Optimization curve of rod lengthρ.Figure 22
Optimal volume.Figure 23
ρ corresponds to the dexterity after optimization.Objective functions are as follows:(33)f1=minρ=2sinβ+cosθcosβ+1L+51,f2=minQD=11+tanθ4tanθ−12.Constraints are as follows:
(34)−π2≤θ≤0;0≤β≤π4;L=100sinβ+153.50≤β≤π4.The parameters of multiobjective optimization algorithm are shown in Table1:Table 1ParetofractionPopulation sizeGenerationsStallGenLimitTolFun0.41002002001e-100As shown in Figure23, ρmin=126.3mm.It is calculated by equations (28), (29), and (31): β=0.502rad, θ=−arcsintanβ=−0.581rad(35)L=100sinβ+153.5=1000.481+153.5=361.7mm.Paretooptimal solutions (8 groups selected) are shown in Table2.Table 2ValuesθradβradρmmQDUnoptimized-1.00000.500000135.0000000.2213Optimized-1.2702610.321750134.2471560.2452-1.5707710.784929127.8282840.5246-1.4607820.391701144.0552990.0015-0.5801300.502304126.2975730.7158-1.1707570.785391133.4965460.2445-0.7707810.785163137.5199300.2001-1.3707830.425692128.7970830.4209-0.9707790.534215130.7971730.3262Based on the results of parameter optimization, and considering the processing and assembly technology, select the fourth group of parameters in Table2. Assembled haptic device is shown in Figure 24.Figure 24
Assembled haptic device.


## 5.1. Paretooptimal Solution of Multiobjective Optimization

The problem of multiobjective optimization can be described as follows:(32)minf1x,f2x,⋯,fmx,s.t.lb≤x≤ub,Aeq∗x=beq,A∗x≤b.fix is the objective function to be optimized, x (lb≤x≤ub): variables to be optimized.Aeq∗x=beq is the linear equality constraint of x, A∗x≤b: linear inequality constraints.lb≤x≤ub is the range of x.


### 5.1.1. Gamultiobj Function

The gamultiobj function is an improved multiobjective optimization algorithm based on NSGA (a variant of NSGA). The emergence of the gamultiobj function provides a good way to solve the multiobjective optimization problem in the MATLAB platform. Organization chart of multiobjective optimization function gamultiobj based on GA is shown in Figure20. In this paper, the length ρ and dexterity QD are taken as the objective functions, and the independent variables are θ,β. The best matching value between the volume space suitable for the optimal dexterity is shown in Figures 21–23.Figure 20
Organization chart of the gamultiobj function.Figure 21
Optimization curve of rod lengthρ.Figure 22
Optimal volume.Figure 23
ρ corresponds to the dexterity after optimization.Objective functions are as follows:(33)f1=minρ=2sinβ+cosθcosβ+1L+51,f2=minQD=11+tanθ4tanθ−12.Constraints are as follows:
(34)−π2≤θ≤0;0≤β≤π4;L=100sinβ+153.50≤β≤π4.The parameters of multiobjective optimization algorithm are shown in Table1:Table 1ParetofractionPopulation sizeGenerationsStallGenLimitTolFun0.41002002001e-100As shown in Figure23, ρmin=126.3mm.It is calculated by equations (28), (29), and (31): β=0.502rad, θ=−arcsintanβ=−0.581rad(35)L=100sinβ+153.5=1000.481+153.5=361.7mm.Paretooptimal solutions (8 groups selected) are shown in Table2.Table 2ValuesθradβradρmmQDUnoptimized-1.00000.500000135.0000000.2213Optimized-1.2702610.321750134.2471560.2452-1.5707710.784929127.8282840.5246-1.4607820.391701144.0552990.0015-0.5801300.502304126.2975730.7158-1.1707570.785391133.4965460.2445-0.7707810.785163137.5199300.2001-1.3707830.425692128.7970830.4209-0.9707790.534215130.7971730.3262Based on the results of parameter optimization, and considering the processing and assembly technology, select the fourth group of parameters in Table2. Assembled haptic device is shown in Figure 24.Figure 24
Assembled haptic device.


## 5.1.1. Gamultiobj Function

The gamultiobj function is an improved multiobjective optimization algorithm based on NSGA (a variant of NSGA). The emergence of the gamultiobj function provides a good way to solve the multiobjective optimization problem in the MATLAB platform. Organization chart of multiobjective optimization function gamultiobj based on GA is shown in Figure20. In this paper, the length ρ and dexterity QD are taken as the objective functions, and the independent variables are θ,β. The best matching value between the volume space suitable for the optimal dexterity is shown in Figures 21–23.Figure 20
Organization chart of the gamultiobj function.Figure 21
Optimization curve of rod lengthρ.Figure 22
Optimal volume.Figure 23
ρ corresponds to the dexterity after optimization.Objective functions are as follows:(33)f1=minρ=2sinβ+cosθcosβ+1L+51,f2=minQD=11+tanθ4tanθ−12.Constraints are as follows:
(34)−π2≤θ≤0;0≤β≤π4;L=100sinβ+153.50≤β≤π4.The parameters of multiobjective optimization algorithm are shown in Table1:Table 1ParetofractionPopulation sizeGenerationsStallGenLimitTolFun0.41002002001e-100As shown in Figure23, ρmin=126.3mm.It is calculated by equations (28), (29), and (31): β=0.502rad, θ=−arcsintanβ=−0.581rad(35)L=100sinβ+153.5=1000.481+153.5=361.7mm.Paretooptimal solutions (8 groups selected) are shown in Table2.Table 2ValuesθradβradρmmQDUnoptimized-1.00000.500000135.0000000.2213Optimized-1.2702610.321750134.2471560.2452-1.5707710.784929127.8282840.5246-1.4607820.391701144.0552990.0015-0.5801300.502304126.2975730.7158-1.1707570.785391133.4965460.2445-0.7707810.785163137.5199300.2001-1.3707830.425692128.7970830.4209-0.9707790.534215130.7971730.3262Based on the results of parameter optimization, and considering the processing and assembly technology, select the fourth group of parameters in Table2. Assembled haptic device is shown in Figure 24.Figure 24
Assembled haptic device.


## 6. Results and Discussion

According to the characteristics of series and parallel mechanisms, a new six-DOF hybrid haptic device was designed combining series mechanism 3PRPaR and parallel mechanism ROBO_003. The DOF of 3PRPaR is 3. In this paper, the freedom form of parallel mechanism ROBO_003 was analyzed based on screw theory. It can be concluded from equations (5) and (6) that the DOF of the moving platform of ROBO_003 is 3, and the DOF form is rotation. That is, moving platform of ROBO_003 rotates around X, Y, and Z axes, respectively. The kinematics of 3PRPaR was analyzed in this paper. The workspace of the haptic device was analyzed in the Monte Carlo method. As shown in Figures 11–14, the workspace of 3PRPaR contains a cube, and it meets the design requirements of the haptic device. The workspace corresponds to the length and volume of the haptic device. In order to make the length and the volume meet the requirements, the optimization indices are selected to optimize the haptic device. Multiobjective optimization algorithm based on GA was selected. θ,β was selected as independent variables, Length of moving rod ρand dexterityQD was selected as objective functions. The optimal angle θ,β is -0.581 rad and 0.502 rad obtained, and L=361.7mm. The optimal length ρ and optimal dexterity QDof the haptic device were obtained. The fourth data is selected from 8 paretooptimal solutions in Table 2. The optimal length is ρmin=126.3mm, and the optimal dexterity is QD=0.7158.Minimum ρ means minimum volume, and minimum volume and optimal dexterity mean optimal performance of the haptic device. The optimized dexterity is 69% higher than unoptimized, and the rod length is reduced by 6%.It can be concluded in Figure 23 that ρand QD are approximately inversely proportional. The larger the volume is, the worse the dexterity of the haptic device.The requirements of the haptic device are tracking performance and force feedback performance, which are determined by the appropriate volume and dexterity. When the rod length is determined, the performance of haptic device is determined by the maximum volume and optimal dexterity of the haptic device. In this paper, the multiobjective optimization algorithm based on the genetic algorithm was used to optimize the relevant parameters to obtain the optimal results. It provides a theoretical basis for the designation of haptic devices.


## 7. Conclusion

A new six-DOF hybrid haptic device was presented in this paper. Screw theory was carried out to analyze the DOF of ROBO_003.The workspace of the haptic device is determined in the Monte Carlo method. Based on the genetic algorithm, the parameters of the mechanism are optimized according to the dexterity index and volume evaluation index, and the results showed that the optimal parameters were obtained. The results meet the designation requirements of the haptic device. As an indispensable part of virtual reality interaction, haptic device has become a hot spot in the research field of robot technology. Haptic devices are applied to teleoperation robots which are widely used in underwater, space, and land dangerous places. In the future, there will be more mathematical theories applied in the theory of mechanism, kinematics, and dynamics. It has great application and research value.



---

*Source: 5551585-2021-09-08.xml*