Performance analysis and multiobjective structural optimization of 4RRR pure rotation parallel mechanism

Abstract

In this article, the performance analysis and multiobjective structure optimization of 4RRR parallel mechanism are carried out. Firstly, the 4RRR pure rotation parallel mechanism and its design route are introduced. Secondly, the Jacobian matrices in 2DoF pure rotation and 3DoF pure rotation modes are derived using the motion equations of the mechanism. Next, the singularity analysis, kinematic dexterity analysis, dynamic dexterity analysis, and stiffness analysis of the mechanism are carried out, respectively, and it is proved that there is no singularity in the mechanism in its workspace. Since the dexterity performance expression is a nonlinear piecewise function, the kinematic local comprehensive dexterity index and the dynamic local comprehensive dexterity index are proposed as the objects of analysis. Furthermore, the kinematic global comprehensive dexterity index, the dynamic global comprehensive dexterity index, and the global comprehensive stiffness index are proposed to carry out the multiobjective structural optimization. Finally, NSGA3 was used to complete the optimization, and the comprehensive optimal solution of the structure size was obtained.

1. Introduction

Parallel mechanism has been widely concerned and studied in industry and academia because of its advantages of high stiffness, high control fault tolerance, easy motion decoupling, and high dexterity [1]. The application of parallel mechanism can be seen in many fields [2, 3]. Especially in solving the problem of complex surface machining, the parallel mechanism shows good adaptability. The serial–parallel hybrid machine tool composed of 2–5DoF parallel mechanism and 1–2DoF tool has become mature at present, and typical products such as Tricept [4], Z3 head [5], Exechon [6, 7], and A3 head [8] have been produced. At the same time, in the medical field, the demand for high stiffness and high precision, as well as the motion requirements of pitch, translation, and rotation, can also be met by parallel mechanisms, such as limb rehabilitation [9] and minimally invasive surgery [10].

However, the problems of small workspace, discontinuous dexterous space, and high disorder of singularity of parallel mechanism are always accompanied by its development process. In order to quantify the performance of parallel mechanism, a series of evaluation indexes are proposed. The kinematic performance analysis of the mechanism can be carried out based on the Jacobian matrix of the mechanism: the condition number index is widely used to evaluate the dexterity and singularity of the parallel mechanism [11–14]. Singularity analysis is very important. When the mechanism is close to the singularity, its kinematic performance and stiffness will be degraded to varying degrees. The controllability index is used to measure the transmission quality of speed and force, which can be divided into the speed controllability index [15, 16] and the force controllability index [17]. The current mainstream direction of dexterity of parallel mechanisms is to analyze kinematic dexterity [11] and dynamic dexterity [18]. The two indexes are the ratio of the minimum eigenvalue to the maximum eigenvalue of the matrix (kinematic dexterity corresponds to the Jacobian matrix, and dynamic dexterity corresponds to the mass matrix) as a reference. Assuming that the minimum eigenvalue is the short axis of the ellipse and the maximum eigenvalue is the long axis of the ellipse, the closer the ellipse is to the positive circle, the higher the dexterity performance is.

At present, the above performance evaluation indexes are widely used in the design and development process of parallel mechanisms in industry and academia [19, 20]. However, it is difficult to achieve the optimal solution of multiple performance indexes at the same time. Therefore, the introduction of multiobjective optimization algorithm to deal with such problems is the current mainstream idea [21]. Genetic algorithm (GA) is an intelligent algorithm that simulates the genetic and evolutionary laws of natural organisms. In the optimization process, the data obtain a new population through evolutionary processes such as repeated selection, crossover, and elite retention until the convergence condition is reached. Yang et al. used GA to optimize the workspace section of the mechanism to improve the working range [22]. Ye et al. used GA to optimize the motion/force transmission index of the mechanism to obtain the optimal size of the mechanism [23]. Wu used GA to optimize stiffness performance, motion performance, workspace, and dynamic performance [24, 25]. Bounab et al. and Kelaiaia et al. optimized the stiffness performance, system mass, and workspace as objective functions based on GA [26,26]. Various derivative algorithms of GA such as NSGA2 and NSGA3 are widely used in academia [27–29]. The particle swarm optimization (PSO) algorithm simulates the foraging behavior of the bird swarm. Each individual in the bird swarm updates its spatial position and flight speed by tracking two extreme values. The optimal particle found by a single individual in the iterative process is the individual extreme value, while the optimal particle found by the group in the iterative process is the global extreme value. Lian used PSO to optimize the stiffness, mass, and dynamic performance of PaQuad parallel mechanism [30]. Qi et al. used PSO algorithm to optimize the workspace, kinematics, stiffness, and dynamic performance of the mechanism considering the size parameters and cross-section parameters of the mechanism [31]. Wang et al. optimized the workspace, input coupling rate, and natural frequency of the mechanism based on PSO [32]. Sun et al. optimized the mass and stiffness of the T5 parallel mechanism based on PSO [33]. Inspired by GA and PSO, various bionic algorithms have emerged, such as multiobjective dragonfly algorithm [34, 35], gazelle optimization algorithm [36, 37], flower pollination algorithm [38, 39], and gray wolf optimization algorithm [40, 41].

In this study, the surface roll press of intricate free-form surfaces serves as the backdrop, with the objective being achieved through a parallel mechanism that carries the end-effector tool. Recognizing the necessity for the mechanism to possess decoupled pure rotation capabilities, a 4RRR parallel mechanism augmented with redundant RRR branches is introduced, derived from the foundational 3RRR mechanism. The tool is mounted at the position of the moving platform within this mechanism, while the static platform is fixedly connected to the mobile unit of the linear motion module. To investigate the impact of the structural dimensions of the 4RRR mechanism on its kinematic performance, it is imperative to conduct analyses on singularity, dexterity, and stiffness. Furthermore, multiobjective optimization is employed to ascertain the optimal mechanism dimensions that yield the most favorable comprehensive performance.

The structure of this article is as follows. In the second section, the optimized 4RRR parallel mechanism is introduced, and the kinematics analysis is completed. The velocity Jacobian matrix is derived according to the kinematics equation of the moving platform. The third section uses the velocity Jacobian matrix to analyze the singularity of the parallel mechanism. In the fourth section, the kinematics dexterity analysis is based on the velocity Jacobian matrix, and the dynamics dexterity analysis is based on the moving platform mass matrix. The fifth section is focuses on stiffness analysis, and the stiffness expression of the parallel mechanism system is completed based on the principle of virtual work. In the sixth section, the multiobjective optimization is carried out with the kinematic global comprehensive dexterity index (KGCDI), the dynamic global comprehensive dexterity index (DGCDI), and global comprehensive stiffness index (GCSI) proposed in this article as the objective function, and the structural size of the parallel mechanism as the optimization variable. NSGA3 is selected as the optimization algorithm. Finally, the seventh section summarizes the full text.

2. Kinematics analysis of 4RRR parallel mechanism

2.1. Research background

This study involves laying down surface substrate on a specialized workpiece and ensuring a tight fit between the surface substrate and the workpiece through roller rolling and pressing. The workpiece exhibits complex local surface morphology features, and the surface of workpiece is shown in Figure 1(a). After several bus bars are taken from the surface of the workpiece and expressed in three-dimensional space, it can be seen that each bus bar is not parallel in the three coordinate planes (as shown in Figure 1(b)), that is, three rotational degrees of freedom are required for movement along the surface of the curved surface.

Figure 1. intricate free-form surfaces of workpiece.

2.2. Design of 4RRR parallel mechanism

The 3RRR mechanism has received extensive attention and research as a pure rotational parallel mechanism. As a 3DoF mechanism, it inherently necessitates three actuators, and theoretically, three single-actuator kinematic chains suffice to achieve its motion capabilities. However, this study innovates in enhancing the load-bearing capacity and motion stability of the 3RRR mechanism by incorporating an additional RRR branched chain as a redundant link (Simulation tests can be seen in Section 5.4). In addition, because the new redundant branch chain has the same RRR structure as the drive branch chain, many parts are shared between the two types of branch chains, and the error of parts in the same batch of processing and assembly is within the same range. Therefore, it can ensure the stability of the mechanism precision of the same batch production. At the same time, it can show good part interchangeability in the process of parts replacement. To ensure that the three rotational degrees of freedom strictly adhere to the mutual perpendicularity of the ‘xyz’ axes, the four kinematic chains must be arranged with 90° intervals.

Figure 2(a)–(c) illustrates the conceptual design approach of the mechanism. Upon obtaining the structural theoretical model depicted in Figure 2(d), it was further materialized into a 3D analytical model shown in Figure 2(e). The four branched chains of the mechanism share an identical structural form, with each chain’s revolute pair axis closest to the moving platform residing within the plane of the moving platform, and adjacent axes being mutually perpendicular. Ultimately, after the structural design and selection of the actuator placement, the prototype model of the 4RRR mechanism was developed, as presented in Figure 2(f). It is worth noting that, unlike the driven chain, the no-drive chain has an electrically controlled locking function to obtain stiffness. When the no-drive branch chain is moved to the target position, the cow-eye wheel will lock at this position.

Figure 2. Structural diagram of rolling parallel mechanism.

2.3. Degree of freedom analysis based on modified GK formula

Furthermore, the modified GK formula (Eq. (1)) is used to analyze the degree of freedom of the parallel mechanism. Because there is a redundant branch chain in the mechanism, it needs to be excluded in the analysis process.

(1) \begin{equation} M=d\left(n-g-1\right)+\sum _{i=1}^{g}f_{i}+v-\zeta \end{equation}

where d is the order of the mechanism (d = 3), n is the number of constructions including the rack (n = 8), g is the number of kinematic pairs (g = 9), f _i is the degree of freed`om of the i ^th kinematic pair ( $\sum _{i=1}^{g}f_{i}$ = 9), v is the redundant degree of freedom (v = 0), and ζ is the local degree of freedom (ζ = 0). Finally, substituting the above parameters into Eq. 1, the degree of freedom of the mechanism can be calculated as M = 3.

2.3. Kinematics analysis of 4RRR parallel mechanism

Since the four chains of the mechanism are essentially identical, the kinematic equations of motion platform can be solved by taking one of the chains. Taking chain 1 as an example, a kinematic diagram of the mechanism shown in Figure 1(e) and a coordinate transformation from point O to point P are established. The {1} is the reference frame, and when the roller at the end point has a feed rate of 0, it is in the {4} frame at point P. When the roller has a feed rate, it is in the {5} frame at point P. The rotational quantities corresponding to the three revolute joints of the chain are denoted as α _i , β _i , and γ _i (where i = 1, 2, 3), with the first revolute joint being the chain driver. The rotations of point P about the x, y, and z axes are denoted as θ ₁, θ ₂, and θ ₃, respectively. Taking chain 1 as an example and utilizing the geometric relationship between the z-axes of frames {2} and {3} being orthogonal, an expression relating the rotational input to the positional output can be derived (Eq. (2))

(2) \begin{align} &\frac{-d_{1}r\times \sin \alpha _{1}\cos \theta _{2}+d_{1}d_{3}\times \sin \theta _{2}}{d_{1}\sqrt{r^{2}+{d_{3}}^{2}}}=0\Rightarrow \sin \alpha _{1}=\frac{d_{3}}{r}\tan \theta _{2}\nonumber\\& \mathit{General}\quad \mathit{Solution}:\quad \sin \alpha _{i}=k\tan \theta _{j}\left(i=1,2;j=2,1;k=\frac{d_{3}}{r}\right)\end{align}

Furthermore, the expression of the relationship between the rotational input and the attitude output of the three driving branches can be obtained. The velocity Jacobian matrix (Eq. (3)) and the acceleration equations of the mechanism (Eq. (4)) can be obtained by differentiating and sorting the two sides of the equation.

(3) \begin{equation} \left[\dot{\alpha }\right]=J_{v}\left[\dot{\theta }\right],\textit{where}\left[\dot{\alpha }\right]=\left[\begin{array}{l} \dot{\alpha }_{1}\\ \dot{\alpha }_{2}\\ \dot{\alpha }_{3} \end{array}\right],\left[\dot{\theta }\right]=\left[\begin{array}{l} \dot{\theta }_{1}\\ \dot{\theta }_{2}\\ \dot{\theta }_{3} \end{array}\right] \end{equation}

(4) \begin{equation} \left[\ddot{\alpha }\right]=\dot{J}_{v}\left[\dot{\theta }\right]+J_{v}\left[\ddot{\theta }\right],\textit{where}\left[\ddot{\alpha }\right]=\left[\begin{array}{l} \ddot{\alpha }_{1}\\ \ddot{\alpha }_{2}\\ \ddot{\alpha }_{3} \end{array}\right],\left[\ddot{\theta }\right]=\left[\begin{array}{l} \ddot{\theta }_{1}\\ \ddot{\theta }_{2}\\ \ddot{\theta }_{3} \end{array}\right] \end{equation}

According to the structural characteristics of the 4RRR mechanism, its Jacobian matrix can be divided into two cases:

I. When the rotation velocity of branched chain 1, branched chain 2, and branched chain 3 are not equal, the mechanism has two degrees of freedom around the x-axis and around the y-axis. At this time, the velocity Jacobian matrix and acceleration equations can be expressed as Eq. (5) and Eq. (6).

(5) \begin{equation} {}^{1}{\boldsymbol{J}}{_{v}^{}}=\left[\begin{array}{ll} \dfrac{\partial \alpha _{1}\left(\theta _{1},\theta _{2}\right)}{\partial \theta _{2}} & \dfrac{\partial \alpha _{1}\left(\theta _{1},\theta _{2}\right)}{\partial \theta _{1}}\\ \dfrac{\partial \alpha _{2}\left(\theta _{1},\theta _{2}\right)}{\partial \theta _{2}} & \dfrac{\partial \alpha _{2}\left(\theta _{1},\theta _{2}\right)}{\partial \theta _{1}} \end{array}\right]=\left[\begin{array}{ll} \dfrac{\dfrac{r}{d_{3}}\cos \theta _{1}}{1+\left(\dfrac{r}{d_{3}}\sin \theta _{1}\right)^{2}} & 0\\ 0 & \dfrac{\dfrac{r}{d_{3}}\cos \theta _{2}}{1+\left(\dfrac{r}{d_{3}}\sin \theta _{2}\right)^{2}} \end{array}\right] \end{equation}

(6) \begin{equation} \left[\ddot{\alpha }\right]=\left[\begin{array}{ll} \dfrac{-\dfrac{r}{d_{3}}\left(3+\dfrac{r}{d_{3}}\right)\sin \theta _{1}\dot{\theta }_{1}}{\left[1+\left(\dfrac{r}{d_{3}}\sin \theta _{1}\right)^{2}\right]^{2}} & 0\\ 0 & \dfrac{-\dfrac{r}{d_{3}}\left(3+\dfrac{r}{d_{3}}\right)\sin \theta _{2}\dot{\theta }_{2}}{\left[1+\left(\dfrac{r}{d_{3}}\sin \theta _{2}\right)^{2}\right]^{2}} \end{array}\right]\left[\dot{\theta }\right]+{}^{1}{\boldsymbol{J}}{_{v}^{}}\left[\ddot{\theta }\right] \end{equation}

II. When the branch 1, branch 2, and branch 3 are at the same rotation velocity, the degree of freedom of the mechanism around the z-axis will be released, and the velocity Jacobian matrix and acceleration equations can be expressed as Eq. (7) and Eq. (8).

(7) \begin{equation} {}^{2}{\boldsymbol{J}}{_{v}^{}}=\left[\begin{array}{lll} \dfrac{\partial \alpha _{1}\left(\theta _{1},\theta _{2},\theta _{3}\right)}{\partial \theta _{2}} & \dfrac{\partial \alpha _{1}\left(\theta _{1},\theta _{2},\theta _{3}\right)}{\partial \theta _{1}} & \dfrac{\partial \alpha _{1}\left(\theta _{1},\theta _{2},\theta _{3}\right)}{\partial \theta _{3}}\\ \dfrac{\partial \alpha _{2}\left(\theta _{1},\theta _{2},\theta _{3}\right)}{\partial \theta _{2}} & \dfrac{\partial \alpha _{2}\left(\theta _{1},\theta _{2},\theta _{3}\right)}{\partial \theta _{1}} & \dfrac{\partial \alpha _{2}\left(\theta _{1},\theta _{2},\theta _{3}\right)}{\partial \theta _{3}}\\ \dfrac{\partial \alpha _{3}\left(\theta _{1},\theta _{2},\theta _{3}\right)}{\partial \theta _{2}} & \dfrac{\partial \alpha _{3}\left(\theta _{1},\theta _{2},\theta _{3}\right)}{\partial \theta _{1}} & \dfrac{\partial \alpha _{3}\left(\theta _{1},\theta _{2},\theta _{3}\right)}{\partial \theta _{3}} \end{array}\right]=\left[\begin{array}{lll} \dfrac{\dfrac{r}{d_{3}}\cos \theta _{1}}{1+\left(\dfrac{r}{d_{3}}\sin \theta _{1}\right)^{2}} & 0 & 0\\ 0 & \dfrac{\dfrac{r}{d_{3}}\cos \theta _{2}}{1+\left(\dfrac{r}{d_{3}}\sin \theta _{2}\right)^{2}} & 0\\ 0 & 0 & 1 \end{array}\right] \end{equation}

(8) \begin{equation} \left[\ddot{\alpha }\right]=\left[\begin{array}{lll} \dfrac{-\dfrac{r}{d_{3}}\left(3+\dfrac{r}{d_{3}}\right)\sin \theta _{1}\dot{\theta }_{1}}{\left[1+\left(\dfrac{r}{d_{3}}\sin \theta _{1}\right)^{2}\right]^{2}} & 0 & 0\\ 0 & \dfrac{-\dfrac{r}{d_{3}}\left(3+\dfrac{r}{d_{3}}\right)\sin \theta _{2}\dot{\theta }_{2}}{\left[1+\left(\dfrac{r}{d_{3}}\sin \theta _{2}\right)^{2}\right]^{2}} & 0\\ 0 & 0 & 1 \end{array}\right]\left[\dot{\theta }\right]+{}^{2}{\boldsymbol{J}}{_{v}^{}}\left[\ddot{\theta }\right] \end{equation}

where θ ₃ is a dummy variable, whose value corresponds to the angle of rotation when θ ₁, and θ ₂ move synchronously.

2.5. Characteristics of the novel mechanism

The 4RRR mechanism adds a redundant branch chain with the same structure but no actuation to the 3RRR mechanism. Increasing the number of branches can improve the overall stiffness of the mechanism. The four branches are equally spaced to distribute the load, especially considering that when there are no redundant branches, each revolute pair of the branch on the opposite side of that position will experience a large torque during the driving process.

3. Singularity analysis of 4RRR parallel mechanism

The discrimination method based on Jacobian matrix is the most commonly used method to solve the singularity of mechanisms. The singularity problems of parallel mechanisms can be divided into three categories: inverse kinematics singularity (IKS), forward kinematics singularity (FKS), and combined singularity (CS). These three kinds of singularities can be distinguished by Jacobian matrix determinant [42–45].

The general formula of the mechanism motion relationship can be expressed as $\dot{\boldsymbol{I}}=\sqrt{\boldsymbol{J}_{s}^{-1}\boldsymbol{J}_{p}}\dot{\boldsymbol{O}}=\boldsymbol{J}\dot{\boldsymbol{O}}$ , where $\dot{\boldsymbol{I}}$ is the input velocity vector matrix, $\dot{\boldsymbol{O}}$ is the output velocity vector matrix, $\boldsymbol{J}=\sqrt{\boldsymbol{J}_{s}^{-1}\boldsymbol{J}_{p}}$ is the mechanism Jacobian matrix, and $\boldsymbol{J}_{s}$ and $\boldsymbol{J}_{p}$ represent the inverse Jacobian matrix and the positive Jacobian matrix, respectively. At this point, the above three mechanism singularity solution methods of IKS, FKS, and CS can be expressed as shown in Eq. (9).

(9) \begin{equation} \left\{\begin{array}{l} \mathrm{IKS}\colon \begin{array}{l} \end{array}\det \left(\boldsymbol{J}_{s}\right)=0,\det \left(\boldsymbol{J}_{ps}\right)\neq 0\\ \mathrm{FKS}\colon \begin{array}{l} \end{array}\det \left(\boldsymbol{J}_{s}\right)\neq 0,\det \left(\boldsymbol{J}_{ps}\right)=0\\ \mathrm{CS}\colon \begin{array}{l} \end{array}\det \left(\boldsymbol{J}_{s}\right)=0,\det \left(\boldsymbol{J}_{ps}\right)=0 \end{array}\right. \end{equation}

$\boldsymbol{J}_{ps}$ is a submatrix of $\boldsymbol{J}_{p}$ , and $\boldsymbol{J}_{ps}=\boldsymbol{J}_{p}$ when $\boldsymbol{J}_{p}$ is a square matrix. The output angle variables of the mechanism are sequentially valued in the reachable space, and the spatial distribution of the determinant of the inverse Jacobian matrix and the positive Jacobian matrix is obtained, as shown in Figure 3, where $\det (\boldsymbol{J}_{s})$ and $\det (\boldsymbol{J}_{p})$ are as follows, respectively.

Figure 3. The determinant of the forward and inverse Jacobian.

(10) \begin{equation} \det \left(\boldsymbol{J}_{s}\right)=\frac{\left[1+\left(\frac{r}{d_{3}}\sin \theta _{1}\right)^{2}\right]\left[1+\left(\frac{r}{d_{3}}\sin \theta _{2}\right)^{2}\right]}{\left(\frac{r}{d_{3}}\right)^{2}\cos \theta _{1}\cos \theta _{2}} \end{equation}

(11) \begin{equation} \det \left(\boldsymbol{J}_{p}\right)=\frac{\left(\frac{r}{d_{3}}\right)^{2}\cos \theta _{1}\cos \theta _{2}}{\left[1+\left(\frac{r}{d_{3}}\sin \theta _{1}\right)^{2}\right]\left[1+\left(\frac{r}{d_{3}}\sin \theta _{2}\right)^{2}\right]} \end{equation}

According to Eq. (5) and Eq. (7), $\det (\boldsymbol{J}_{s})$ is always greater than 0, so there is no IKS space and CS space. When $\theta _{1}=\theta _{2}=\pi /2$ , there is a unique point in the FKS space. However, the reachable workspace of the mechanism rotating around two axes is $-26.6^{\circ}\leq \theta _{1}\boldsymbol{,}\theta _{2}\leq 26.6^{\circ}$ , and there are no singular points in this range. Therefore, there is no singularity in the workspace of the mechanism.

4. Dexterity analysis of 4RRR parallel mechanism

4.1. Kinematics dexterity analysis

The kinematics dexterity of parallel mechanism is to measure the ability to perform various actions in its workspace, which is usually measured from the workspace range, motion degree of freedom, and obstacle avoidance ability. Salisbury et al. proposed to use the condition number of the Jacobian matrix to evaluate the kinematics dexterity of the mechanism [46].

The condition number of the matrix is the product of the matrix spectral norm and the inverse matrix spectral norm: $\text{cond}(\boldsymbol{A})=\kappa (\boldsymbol{A})=\| \boldsymbol{A}\| \cdot \| \boldsymbol{A}^{-1}\|$ , and the condition number of the Jacobian matrix can be expressed as Eq. (12).

(12) \begin{equation} \text{cond}\left(\boldsymbol{J}\right)=\kappa \left(\boldsymbol{J}\right)=\left\| \boldsymbol{J}\right\| \cdot \left\| \boldsymbol{J}^{-1}\right\| \end{equation}

The spectral norm is defined as the maximum singular value of the matrix, so there is the following relationship.

(13) \begin{equation} \left\{\begin{array}{l} \left\| \boldsymbol{J}\right\| =\sigma _{\boldsymbol{\max }}\\ \left\| \boldsymbol{J}^{-1}\right\| =\frac{1}{\sigma _{\boldsymbol{\min }}} \end{array}\right. \end{equation}

(14) \begin{equation} 1\leq \kappa \left(\boldsymbol{J}\right)=\left\| \boldsymbol{J}\right\| \cdot \left\| \boldsymbol{J}^{-1}\right\| =\frac{\sigma _{\boldsymbol{\max }}}{\sigma _{\boldsymbol{\min }}}\lt \infty \end{equation}

When evaluating the kinematics local dexterity, $1\leq \kappa (\boldsymbol{J})\lt \infty$ or $0\lt 1/\kappa (\boldsymbol{J})\leq 1$ can be used as the evaluation index. The closer $\kappa (\boldsymbol{J})$ is to 1 means that the dexterity of the mechanism is higher, and the closer ∞ means that the position of the mechanism is closer to singularity. On the contrary, $1/\kappa (\boldsymbol{J})$ close to 1 indicates that the higher the dexterity of the mechanism, the closer to 0 indicates that the position of the mechanism is closer to the singularity. In this article, $1/\kappa (\boldsymbol{J})$ is used as the dexterity evaluation index. By setting the structural parameter $k=r/d_{3}$ , the expression of kinematics local dexterity (Eq. (15)) can be derived.

(15) \begin{equation} 1/\kappa \left(\boldsymbol{J}\right)=\dfrac{\sigma _{\boldsymbol{\min }}}{\sigma _{\boldsymbol{\max }}}=\left\{\begin{array}{l} \dfrac{k\cos \theta _{1}}{1+\left(k\sin \theta _{1}\right)^{2}}\boldsymbol{/}\dfrac{k\cos \theta _{2}}{1+\left(k\sin \theta _{2}\right)^{2}}\begin{array}{l} \end{array}\begin{array}{l} \end{array}\begin{array}{l} \end{array}s\mathit{.}t\mathit{.}\dfrac{k\cos \theta _{1}}{1+\left(k\sin \theta _{1}\right)^{2}}\leq \dfrac{k\cos \theta _{2}}{1+\left(k\sin \theta _{2}\right)^{2}}\\ \dfrac{k\cos \theta _{2}}{1+\left(k\sin \theta _{2}\right)^{2}}\boldsymbol{/}\dfrac{k\cos \theta _{1}}{1+\left(k\sin \theta _{1}\right)^{2}}\begin{array}{l} \end{array}\begin{array}{l} \end{array}\begin{array}{l} \end{array}s\mathit{.}t\mathit{.}\dfrac{k\cos \theta _{1}}{1+\left(k\sin \theta _{1}\right)^{2}}\gt \dfrac{k\cos \theta _{2}}{1+\left(k\sin \theta _{2}\right)^{2}} \end{array}\right. \end{equation}

When k > 1, the shape of the mechanism is flat. When k < 1, the shape of the mechanism is thin and tall. When k = 1, the rotation radius of the first axis of the branch chain of the mechanism is equal to the height of the moving platform. Taking k = 4,2,1 and 0.5 as examples, the distribution of kinematics dexterity in the workspace is shown in Figure 4. The smaller the value of k, the wider the distribution of dexterous space in the workspace.

Figure 4. Kinematics dexterity space of the mechanism under different values of k.

4.2. Dynamics dexterity analysis

The dynamics dexterity is used to measure the acceleration performance of the mechanism at a specific position. Cui (2015) et al. established the dynamics dexterity index using the generalized inertial ellipse (GIE) proposed by Asada (1983) to determine the dynamics dexterity [12, 20]. The mass matrix M _n×n of the mechanism is used to represent the GIE in n-dimensional space. The long axis of the ellipse is the square root of the maximum eigenvalue of M , which is $\sqrt{\lambda _{\boldsymbol{\max }}}$ . The short axis is the square root of the minimum eigenvalue of M , which is $\sqrt{\lambda _{\boldsymbol{\min }}}$ . When the ellipsoid is closer to the sphere, the dynamics performance of the mechanism is better. Therefore, the dynamics local dexterity index (LDI) can be expressed as Eq. (16), and the closer $1/\kappa (\boldsymbol{M})$ is to 1, the better the dynamics dexterity.

(16) \begin{equation} 0\lt 1/\kappa \left(\boldsymbol{M}\right)=\frac{\lambda _{\boldsymbol{\min }}}{\lambda _{\boldsymbol{\max }}}\leq 1 \end{equation}

The mass matrix of the mechanism can be obtained from the dynamic equation of the end of the mechanism. The Lagrange method is to establish the dynamic equation from system energy, and it is necessary to solve the end kinetic energy and potential energy. The mechanism is a 3DoF pure rotation parallel mechanism, so the moving platform motion only has angular velocity but not linear velocity. The structure and size of the moving platform are shown in Figure 5. The moment of inertia of the moving platform can be obtained as shown in Eqs. (17)-(18).

Figure 5. Moving platform model structure and size.

(17) \begin{equation} \left\{\begin{array}{l} I_{\unicode{x2460}}=\frac{1}{2}m_{1}r_{1}^{2}=\frac{\pi }{2}r_{1}^{4}l_{1}\rho _{1}\\ I_{\unicode{x2461}}=m_{2-1}a_{1}^{2}-m_{2-2}a_{2}^{2}=\left(a_{1}^{4}-a_{2}^{4}\right)h\rho _{1}\\ I_{\unicode{x2462}}=\frac{1}{2}m_{3}r_{1}^{2}=\frac{\pi }{2}r_{1}^{4}l_{3}\rho _{2}\\ I_{\unicode{x2463}}=\frac{1}{2}m_{4}r_{2}^{2}=\frac{\pi }{2}r_{2}^{4}l_{2}\rho _{3} \end{array}\right. \end{equation}

(18) \begin{equation} \begin{array}{ll} I & =4\times I_{\unicode{x2460}}+I_{\unicode{x2461}}+2\times I_{\unicode{x2462}}+I_{\unicode{x2463}}\\ & =2\pi r_{1}^{4}l_{1}\rho _{1}+\left(a_{1}^{4}-a_{2}^{4}\right)h\rho _{1}+\pi r_{1}^{4}l_{3}\rho _{2}+\frac{\pi }{2}r_{2}^{4}l_{2}\rho _{3} \end{array} \end{equation}

where $\rho_1$ is the density of $\unicode{x2460}$ and $\unicode{x2461}$ , $\rho_2$ is the density of $\unicode{x2462}$ , and $\rho_3$ is the density of $\unicode{x2463}$ . Further, the kinetic energy of the moving platform can be derived as shown in Eqs. (19)-(20).

(19) \begin{equation} \begin{array}{c} E_{k}=\sum _{i=1}^{3}\frac{1}{2}I\dot{\boldsymbol{\theta }}_{i}^{2}=\sum _{i=1}^{3}\frac{1}{2}I\boldsymbol{J}_{i\boldsymbol{.}}^{-1}\dot{\boldsymbol{\alpha }}^{2}\\ \end{array} \end{equation}

(20) \begin{equation} E_{k}=\frac{1}{2}\left[2\pi r_{1}^{4}l_{1}\rho _{1}+\left(a_{1}^{4}-a_{2}^{4}\right)h\rho _{1}+\pi r_{1}^{4}l_{3}\rho _{2}+\frac{\pi }{2}r_{2}^{4}l_{2}\rho _{3}\right]\left[\frac{1+\left(k\sin \theta _{1}\right)^{2}}{k\cos \theta _{1}}{\dot{\alpha }_{1}}^{2}+\frac{1+\left(k\sin \theta _{2}\right)^{2}}{k\cos \theta _{2}}{\dot{\alpha }_{2}}^{2}+{\dot{\alpha }_{3}}^{2}\right] \end{equation}

The moving platform is only subjected to gravitational potential energy, so the potential energy can be expressed as Eqs. (21)–(22).

(21) \begin{equation} E_{p}=\left(4m_{1}+m_{2}+2m_{3}+m_{4}\right)gd_{4} \end{equation}

(22) \begin{equation} E_{p}=\left[4\pi r_{1}^{2}l_{1}\rho _{1}+\left(a_{1}^{2}-a_{2}^{2}\right)h\rho _{1}+2\pi r_{1}^{2}l_{3}\rho _{2}+\pi r_{2}^{2}l_{2}\rho _{3}\right]gd_{4} \end{equation}

Then, the Lagrange function can be derived as shown in Eqs. (23)–(24).

(23) \begin{equation} L=E_{k}-E_{p} \end{equation}

(24) \begin{equation} \begin{array}{c} L=\frac{1}{2}\left[2\pi r_{1}^{4}l_{1}\rho _{1}+\left(a_{1}^{4}-a_{2}^{4}\right)h\rho _{1}+\pi r_{1}^{4}l_{3}\rho _{2}+\frac{\pi }{2}r_{2}^{4}l_{2}\rho _{3}\right]\left[\dfrac{1+\left(k\sin \theta _{1}\right)^{2}}{k\cos \theta _{1}}{\dot{\alpha }_{1}}^{2}+\dfrac{1+\left(k\sin \theta _{2}\right)^{2}}{k\cos \theta _{2}}{\dot{\alpha }_{2}}^{2}+{\dot{\alpha }_{3}}^{2}\right]\\ -\left[4\pi r_{1}^{2}l_{1}\rho _{1}+\left(a_{1}^{2}-a_{2}^{2}\right)h\rho _{1}+2\pi r_{1}^{2}l_{3}\rho _{2}+\pi r_{2}^{2}l_{2}\rho _{3}\right]gd_{4} \end{array} \end{equation}

Finally, according to Eq. (25), the dynamic equation of the end of the mechanism is solved and sorted out as shown in Eq. (26). Observing the structure of the dynamic equation, it can be seen that the end system of the mechanism is only affected by the inertial force.

(25) \begin{equation} F_{i}=\frac{d}{dt}\frac{\partial L}{\partial \dot{\alpha }_{i}}-\frac{\partial L}{\partial \alpha _{i}}\begin{array}{c} \end{array}\left(i=1,2,3\right) \end{equation}

(26) \begin{equation} F=\left[\begin{array}{l} F_{1}\\ F_{2}\\ F_{3} \end{array}\right]=\left[\begin{array}{ccc} I\left[1+\left(k\sin \theta _{1}\right)^{2}\right]/\left(k\cos \theta _{1}\right) & 0 & 0\\ 0 & I\left[1+\left(k\sin \theta _{2}\right)^{2}\right]/\left(k\cos \theta _{2}\right) & 0\\ 0 & 0 & I \end{array}\right]\left[\begin{array}{l} \ddot{\alpha }_{1}\\ \ddot{\alpha }_{2}\\ \ddot{\alpha }_{3} \end{array}\right] \end{equation}

The mass matrix M is obtained, and then, the dynamics dexterity condition number can be solved as shown in Eq. (27).

(27) \begin{equation} \boldsymbol{M}=\left[\begin{array}{ccc} I\left[1+\left(k\sin \theta _{1}\right)^{2}\right]/\left(k\cos \theta _{1}\right) & 0 & 0\\ 0 & I\left[1+\left(k\sin \theta _{2}\right)^{2}\right]/\left(k\cos \theta _{2}\right) & 0\\ 0 & 0 & I \end{array}\right] \end{equation}

Taking different k values for example analysis, it can be seen that when k≥2, the best position of dynamics dexterity performance is located in the four corners of the rectangular workspace. At this time, the larger the value of k is, the worse the overall dynamics dexterity performance of the mechanism is. When 2 > k > 1, the best position of dynamics dexterity performance gradually moves from the four corners to the center with the decrease of k value, and the area of the best performance area gradually increases. In this range of values, the overall dynamics dexterity performance of the mechanism is excellent. When k ≤ 1, the best position of the dynamics dexterity performance of the mechanism is located at the center of the rectangular workspace, and the overall performance of the mechanism is optimal when k = 1. With the decrease of k value, the overall performance of the mechanism is degraded.

(28) \begin{equation} \frac{1}{\kappa \left(\boldsymbol{M}\right)}=\left\{\begin{array}{l} \left(k\cos \theta _{1}\right)/\left[1+\left(k\sin \theta _{1}\right)^{2}\right]\\[3pt][3pt] \begin{array}{l} \begin{array}{l} \end{array} \end{array}\begin{array}{l} \end{array}\begin{array}{l} \end{array}s\boldsymbol{.}t\boldsymbol{.}\left[1+\left(k\sin \theta _{1}\right)^{2}\right]/k\cos \theta _{1}\geq \left[1+\left(k\sin \theta _{2}\right)^{2}\right]/k\cos \theta _{2}\geq 1\\[3pt] \left(k\cos \theta _{2}\right)/\left[1+\left(k\sin \theta _{2}\right)^{2}\right]\\[3pt] \begin{array}{l} \end{array}\begin{array}{l} \end{array}\begin{array}{l} \end{array}s\boldsymbol{.}t\boldsymbol{.}\left[1+\left(k\sin \theta _{2}\right)^{2}\right]/k\cos \theta _{2}\geq \left[1+\left(k\sin \theta _{1}\right)^{2}\right]/k\cos \theta _{1}\geq 1\\[3pt] \left[1+\left(k\sin \theta _{2}\right)^{2}\right]\cos \theta _{1}/\left[1+\left(k\sin \theta _{1}\right)^{2}\right]\cos \theta _{2}\\[3pt] \begin{array}{l} \end{array}\begin{array}{l} \end{array}\begin{array}{l} \end{array}s\boldsymbol{.}t\boldsymbol{.}\left[1+\left(k\sin \theta _{1}\right)^{2}\right]/k\cos \theta _{1}\geq 1\geq \left[1+\left(k\sin \theta _{2}\right)^{2}\right]/k\cos \theta _{2}\\[3pt] \left[1+\left(k\sin \theta _{1}\right)^{2}\right]\cos \theta _{2}/\left[1+\left(k\sin \theta _{2}\right)^{2}\right]\cos \theta _{1}\\[3pt] \begin{array}{l} \end{array}\begin{array}{l} \end{array}\begin{array}{l} \end{array}s\boldsymbol{.}t\boldsymbol{.}\left[1+\left(k\sin \theta _{2}\right)^{2}\right]/k\cos \theta _{2}\geq 1\geq \left[1+\left(k\sin \theta _{1}\right)^{2}\right]/k\cos \theta _{1}\\[3pt] \left[1+\left(k\sin \theta _{1}\right)^{2}\right]/\left(k\cos \theta _{1}\right)\\[3pt] \begin{array}{l} \begin{array}{l} \end{array}\begin{array}{l} \end{array}\begin{array}{l} \end{array}s\boldsymbol{.}t\boldsymbol{.}1\geq \left[1+\left(k\sin \theta _{2}\right)^{2}\right]/k\cos \theta _{2}\geq \left[1+\left(k\sin \theta _{1}\right)^{2}\right]/k\cos \theta _{1} \end{array}\\[3pt] \left[1+\left(k\sin \theta _{2}\right)^{2}\right]/\left(k\cos \theta _{2}\right)\\[3pt] \begin{array}{l} \begin{array}{l} \end{array}\begin{array}{l} \end{array}\begin{array}{l} \end{array}s\boldsymbol{.}t\boldsymbol{.} \end{array}1\geq \left[1+\left(k\sin \theta _{1}\right)^{2}\right]/k\cos \theta _{1}\geq \left[1+\left(k\sin \theta _{2}\right)^{2}\right]/k\cos \theta _{2} \end{array}\right. \end{equation}

5. Stiffness analysis of 4RRR parallel mechanism

5.1. Input–output force/torque mapping relationship

Figure 6. Dynamics dexterity space of mechanism under different values of k.

Firstly, the mapping relationship between driving force/torque and output force/torque is established. Under the assumption that the stiffness of the structural part is much larger than the transmission stiffness and the joint friction is not considered, the dynamic equation can be established by the virtual work principle to obtain the input–output mapping relationship modeling. Let $\boldsymbol{F}=[f_{1},f_{2},\ldots, f_{n}]^{T}$ denote the output force / torque vector of the moving platform, $\delta \boldsymbol{x}=[\delta x_{1},\delta x_{2},\ldots, \delta x_{n}]^{T}$ denote the virtual displacement vector of the moving platform, $\boldsymbol{\tau }=[\tau _{1},\tau _{2},\ldots, \tau _{m}]^{T}$ denote the driving force / torque vector, and $\delta \boldsymbol{q}=[\delta q_{1},\delta q_{2},\ldots, \delta q_{m}]^{T}$ denote the virtual displacement vector of the driving joint. At this time, the work done by the driving force/torque is all converted into the work output by the moving platform:

(29) \begin{equation} \boldsymbol{\tau }^{T}\delta \boldsymbol{q}=\boldsymbol{F}^{T}\delta \boldsymbol{x} \end{equation}

According to the kinematics analysis, there is a mapping relationship between the input virtual displacement and the output virtual displacement, which is related to the velocity Jacobian matrix:

(30) \begin{equation} \delta \boldsymbol{q}=\boldsymbol{J}_{v}\delta \boldsymbol{x} \end{equation}

The following equality relations can be obtained by combining Eq. (29) and Eq. (30):

(31) \begin{equation} \frac{\boldsymbol{F}^{T}}{\boldsymbol{\tau }^{T}}=\boldsymbol{J}_{v} \end{equation}

Further derivation can be obtained:

(32) \begin{equation} \boldsymbol{F}={\boldsymbol{J}_{v}}^{T}\boldsymbol{\tau } \end{equation}

It can be seen from the above Eqs that the mapping variable between the output force/torque and the input force/torque of the moving platform is the transpose of the velocity Jacobian matrix.

5.2. Establishment of stiffness model

The stiffness matrix describes the stiffness performance of the system well. The spring stiffness generated by the self-locking of the driving unit and the external load stiffness of the structural part are the two most significant components of the stiffness of the mechanism. $\boldsymbol{N}_{i}=[N_{i1},N_{i2},\ldots, N_{in}]^{T}$ is defined as the external force / torque of the i ^th branch chain of the parallel mechanism, $\Delta \boldsymbol{q}_{i}=[\Delta q_{i1},\Delta q_{i2},\ldots, \Delta q_{in}]^{T}$ as the corresponding elastic deformation of the i ^th branch chain under the action of external force / torque, $k_{\chi i}$ as the spring stiffness of the i ^th branch chain, and $k_{Ei}$ as the external load stiffness of the i ^th branch chain. Therefore, the corresponding relationship between the external force / torque of the branch chain and the elastic deformation of the branch chain can be expressed as

(33) \begin{equation} \boldsymbol{N}_{i}=\left(k_{\chi i}+k_{\varepsilon i}\right)\Delta \boldsymbol{q}_{i} \end{equation}

Continuing to define $\boldsymbol{N}=[\boldsymbol{N}_{1},\boldsymbol{N}_{2},\ldots, \boldsymbol{N}_{n}]^{T}$ as the force/torque of all branches of the parallel mechanism, $\Delta \boldsymbol{q}=[\Delta \boldsymbol{q}_{1},\Delta \boldsymbol{q}_{2},\ldots, \Delta \boldsymbol{q}_{n}]^{T}$ as the corresponding elastic deformation of each branch of the parallel mechanism under the action of external force / torque. The systematic expression of Eq. (33) can be obtained as follows:

(34) \begin{equation} \boldsymbol{N}=\left(\boldsymbol{\chi }+\boldsymbol{\varepsilon }\right)\Delta \boldsymbol{q} \end{equation}

In Eq. (34), $\boldsymbol{\chi }=\text{diag}(\chi _{1},\chi _{2},\ldots, \chi _{m})$ is the spring stiffness of the mechanis $\boldsymbol{\varepsilon }=\text{diag}(\varepsilon _{1},\varepsilon _{2},\ldots, \varepsilon _{m})$ is the external load stiffness of the mechanism.

From Eq. (30), it can be seen that there is a relationship between the output deformation of the moving platform of the parallel mechanism and the input deformation of the driver:

(35) \begin{equation} \Delta \boldsymbol{q}=\boldsymbol{J}_{v}\Delta \boldsymbol{x} \end{equation}

It can be seen from Eq. (32) that the corresponding relationship between the output force/torque and the input force/torque of the parallel mechanism is

(36) \begin{equation} \boldsymbol{F}={\boldsymbol{J}_{v}}^{T}\boldsymbol{N} \end{equation}

Simultaneous Eqs. (34)-(36), we can obtain

(37) \begin{equation} \boldsymbol{F}={\boldsymbol{J}_{v}}^{T}\left(\boldsymbol{\chi }+\boldsymbol{\varepsilon }\right)\boldsymbol{J}_{v}\Delta \boldsymbol{x} \end{equation}

$\boldsymbol{K}_{M}={\boldsymbol{J}_{v}}^{T}(\boldsymbol{\chi }+\boldsymbol{\varepsilon })\boldsymbol{J}_{v}$ is defined as the stiffness matrix of parallel mechanism in generalized coordinates. The stiffness matrix $\boldsymbol{K}_{M}$ of the mechanism studied in this article satisfies the characteristic $\boldsymbol{K}_{M}=\text{diag}(K_{11},K_{22},\ldots, K_{mm})$ , so the stiffness of the parallel mechanism system can be expressed as

(38) \begin{equation} K=\left({K_{11}}^{-1}+{K_{22}}^{-1}+\ldots +{K_{mm}}^{-1}\right)^{-1} \end{equation}

5.3. System stiffness expression

The spring stiffness of the drive motor is one of the important factors to ensure that the mechanism maintains the current pose, so this part needs to be considered when solving the system stiffness. According to the introduction of the selected motor technical manual and the actual measurement and calculation results, the spring stiffness of the drive motor can be set to 10⁵N/m.

The external load stiffness $\varepsilon (k_{F},k_{M})$ includes tension and compression stiffness $k_{F}$ and torsion stiffness $k_{M}$ . As shown in Figure 7, the slipway in the branch chain can be regarded as a rigid body, and the connecting rod can be regarded as an L-shaped elastomer composed of L ₁ and L ₂ segments. A single-branched chain is subjected to a quarter of the external load (G) of the moving platform. At this time, the L ₁ segment is subjected to compressive stress and bending stress, and the L ₂ segment is subjected to bending stress.

Figure 7. Branched chain structure size.

The tensile and compressive stiffness of the branched link comes from the tension/compression of the L ₁ and L ₂ segments, so $k_{F}$ can be expressed as follows according to Hooke’s law:

(39) \begin{equation} k_{F}=\frac{G\cos \theta _{1/2}\cos \psi _{2}}{4\Delta L_{1}}+\frac{G\sin \theta _{1/2}}{4\Delta L_{2}} \end{equation}

(40) \begin{equation} \left\{\begin{array}{ll} \Delta L_{1} & =\dfrac{GL_{1,1}\cos \theta _{1/2}\cos \psi _{2}}{4E_{ABS}D_{1}D_{2}}+\dfrac{GL_{1,2}\cos \theta _{1/2}\cos \psi _{2}}{4E_{ABS}\left(D_{1}D_{2}+2D_{3}D_{4}\right)}\\[10pt] \Delta L_{2} & =\dfrac{GL_{2}\sin \theta _{1/2}}{4E_{ABS}\left(D_{1}D_{2}+2D_{3}D_{4}\right)} \end{array}\right. \end{equation}

In above equations, $\Delta L_{1}$ is the deformation of $L_{1}$ under compression load, $\Delta L_{2}$ is the deformation of $L_{2}$ under tensile load, $L_{1,1}$ is the length of the rectangular section, $L_{1,2}$ is the length of the ‘+’ section, and $\theta _{1/2}$ is the output angle of the moving platform around the x-axis or y-axis.

The expressions obtained by the simultaneous calculation of Eqs. (39) and (40):

(41) \begin{equation} k_{F}=\left[\frac{L_{1,1}}{E_{ABS}D_{1}D_{2}}+\frac{L_{1,2}}{E_{ABS}\left(D_{1}D_{2}+2D_{3}D_{4}\right)}\right]^{-1}+\frac{E_{ABS}\left(D_{1}D_{2}+2D_{3}D_{4}\right)}{L_{2}} \end{equation}

The bending stiffness of the branch chain connecting rod comes from the bending moment of the segment and the segment. In order to simplify the calculation, the bending cross sections are regarded as ‘+’ cross sections, which can be further expressed as

(42) \begin{equation} k_{M}=\frac{{}^{1}{M}{_{L1}^{}}}{\Delta \gamma _{1}}+\frac{{}^{2}{M}{_{L1}^{}}}{\Delta \gamma _{2}}+\frac{M_{L2}}{\Delta \gamma _{3}}=\frac{GL_{1}\cos \theta _{1/2}\cos \psi _{1}}{4\Delta \gamma _{1}}+\frac{GL_{1}\sin \theta _{1/2}}{4\Delta \gamma _{2}}+\frac{GL_{2}\cos \theta _{1/2}\cos \psi _{2}}{4\Delta \gamma _{3}} \end{equation}

(43) \begin{equation} \left\{\begin{array}{ll} \Delta \gamma _{1} & =\dfrac{G{L_{1}}^{2}\cos \theta _{1/2}\cos \psi _{1}}{8E_{ABS}\left(\dfrac{{D_{1}}^{3}D_{2}}{12}+\dfrac{D_{3}{D_{4}}^{3}}{6}+\left(D_{2}+D_{3}\right)D_{3}D_{4}\right)}\\ \Delta \gamma _{2} & =\dfrac{G{L_{1}}^{2}\sin \theta _{1/2}}{8E_{ABS}\left(\dfrac{D_{1}{D_{2}}^{3}}{12}+\dfrac{{D_{3}}^{3}D_{4}}{6}+\left(D_{2}+D_{3}\right)D_{3}D_{4}\right)}\\ \Delta \gamma _{3} & =\dfrac{G{L_{2}}^{2}\cos \theta _{1/2}\cos \psi _{2}}{8E_{ABS}\left(\dfrac{{D_{1}}^{3}D_{2}}{12}+\dfrac{D_{3}{D_{4}}^{3}}{6}+\left(D_{2}+D_{3}\right)D_{3}D_{4}\right)} \end{array}\right. \end{equation}

In the formula, $\Delta \gamma _{1}$ is the cross-section rotation angle of segment $L_{1}$ under bending moment, ${}^{1}{M}{_{L1}^{}}, \Delta \gamma _{2}$ is the cross-section rotation angle of segment $L_{1}$ under bending moment, ${}^{2}{M}{_{L1}^{}}, \Delta \gamma _{3}$ is the cross-section rotation angle of segment $L_{2}$ under bending moment, $M_{L2}, \frac{{D_{1}}^{3}D_{2}}{12}$ is the moment of inertia of the large rectangular area in the middle of the ‘+’ section, and $\frac{{D_{3}}^{3}D_{4}}{6}+\left(D_{2}+D_{3}\right)D_{3}D_{4}$ is the moment of inertia of the small rectangular areas on both sides of the ‘+’ section.

The expression of $k_{M}$ can be obtained by simultaneous calculation of Eqs. (40) and (41):

(44) \begin{equation} \begin{array}{l} k_{M}=\dfrac{2E_{ABS}\left(\dfrac{D_{1}^{3}D_{2}}{12}+\dfrac{D_{3}D_{4}^{3}}{6}+\left(D_{2}+D_{3}\right)D_{3}D_{4}\right)}{L_{1}}+\dfrac{2E_{ABS}\left(\dfrac{D_{1}D_{2}^{3}}{12}+\dfrac{D_{3}^{3}D_{4}}{6}+\left(D_{2}+D_{3}\right)D_{3}D_{4}\right)}{L_{1}}\\ \begin{array}{ll} & \end{array}+\dfrac{2E_{ABS}\left(\dfrac{D_{1}^{3}D_{2}}{12}+\dfrac{D_{3}D_{4}^{3}}{6}+\left(D_{2}+D_{3}\right)D_{3}D_{4}\right)}{L_{2}} \end{array} \end{equation}

For the parallel mechanism with four branches, three drive branches have both spring stiffness and external load stiffness, and the redundant branches without drive only have external load stiffness. Therefore, the expressions of $\boldsymbol{\chi }=\text{diag}(\chi _{1},\chi _{2},\ldots, \chi _{m})$ and $\boldsymbol{\varepsilon }=\text{diag}(\varepsilon _{1},\varepsilon _{2},\ldots, \varepsilon _{m})$ can be obtained:

(45) \begin{equation} \boldsymbol{\chi }=\left[\begin{array}{llll} k_{\chi 1} & & & \\ & k_{\chi 2} & & \\ & & k_{\chi 3} & \\ & & & 0 \end{array}\right] \end{equation}

(46) \begin{equation} \boldsymbol{\varepsilon }=\left[\begin{array}{llll} k_{F1}+k_{M1} & & & \\ & k_{F2}+k_{M2} & & \\ & & k_{F3}+k_{M3} & \\ & & & k_{F4}+k_{M4} \end{array}\right] \end{equation}

For the 3DoF decoupled pure rotational parallel mechanism, the velocity Jacobian matrix is third order, so the Jacobian matrix needs to be augmented to fourth order as shown in Eq. (47) when establishing the stiffness model.

(47) \begin{equation} \boldsymbol{J}_{aug.}=\left[\begin{array}{ll} \boldsymbol{J}_{v3\times 3} & 0\\ 0 & 1 \end{array}\right]=\left[\begin{array}{llll} \dfrac{k\cos \theta _{1}}{1+\left(k\sin \theta _{1}\right)^{2}} & & & \\ & \dfrac{k\cos \theta _{2}}{1+\left(k\sin \theta _{2}\right)^{2}} & & \\ & & 1 & \\ & & & 1 \end{array}\right] \end{equation}

At this time, the stiffness matrix under the generalized coordinates of the mechanism can be expressed as

(48) \begin{equation} \begin{array}{ll} \boldsymbol{K}_{M} & ={\boldsymbol{J}_{aug.}}^{T}\left(\boldsymbol{\chi }+\boldsymbol{\varepsilon }\right)\boldsymbol{J}_{aug.}\\ & =\left[\begin{array}{llll} \boldsymbol{J}_{aug.}\left(1,1\right)^{2}\left(k_{\chi 1}+k_{F1}+k_{M1}\right) & & & \\ & \boldsymbol{J}_{aug.}\left(2,2\right)^{2}\left(k_{\chi 2}+k_{F2}+k_{M2}\right) & & \\ & & \left(k_{\chi 3}+k_{F3}+k_{M3}\right) & \\ & & & k_{F4}+k_{M4} \end{array}\right] \end{array} \end{equation}

To investigate the influence of structural dimension k on the stiffness performance of parallel mechanisms, various values of k were assigned as 0.5, 1, 1.2, 1.5, 2, and 4, respectively, yielding the global stiffness distribution as depicted in Figure 8. It was observed that an increase in the value of structural dimension k led to a higher stiffness value under the motionless state; however, this was accompanied by a widening gap between the maximum and minimum values of the global stiffness. Notably, when k was set to 4, although the maximum stiffness among all tested groups was achieved under the motionless condition, the minimum stiffness value fell below that corresponding to k = 1.5 at the same position.

Figure 8. Global stiffness distribution under different k values.

5.4. Simulation loading analysis

In order to compare the stiffness performance of the 4RRR and 3RRR mechanism, the finite element method is used to explore the influence of mechanism branch chain layout on stability. Figure 9(a), (b) (c), respectively, shows the branch chain stress distribution of the above three mechanisms at the moving platform position under the external load of 10 N along the negative direction of z-axis, 7.07 N along the negative direction of z-axis superimposed 7.07 N along the positive direction of x-axis (The net force is 10 N) and 7.07 N along the negative direction of z-axis superimposed 7.07 N along the positive direction of y-axis (The net force is 10N). The colored area is the stress value range of each branch chain of a single mechanism, and the colored solid line is the average value of the branch chain stress. It can be seen that the stress value of each branch chain of the mechanism with equally spaced branches is the same when subjected to vertical load. When subjected to bias load, the three mechanisms all showed uneven loading of the branches. Among them, the 4RRR mechanism has the strongest stability, the 3RRR mechanism with equally spaced branches is the second, and the 3RRR mechanism with decoupled rotation has the worst stability.

Figure 9. Finite element loading simulation test.

Therefore, the following conclusion can be further obtained: if we want to make 3RRR have decoupling pure rotation ability, we can not only change the distribution spacing of branch chains but also need to increase redundant branch chains to obtain the best performance.

6. Multiobjective structure optimization

Based on the analysis of kinematic and dynamic dexterity, it can be observed that a decrease in the structural parameter k leads to a monotonic increase in kinematic dexterity, whereas the dynamic dexterity performance exhibits a maximum point. Furthermore, the stiffness performance does not necessarily improve with the increase in the value of k. Consequently, it is necessary to conduct a multiobjective optimization that incorporates the aforementioned three performance indices to obtain the optimal structural parameters of the mechanism, ensuring optimal overall performance.

6.1. Determination of objective function, optimization variables,and variable boundary conditions

6.1.1. Proposal of kinematics global comprehensive dexterity index

In the multiobjective optimization of dexterity, the global dexterity index (GDI) is usually calculated according to the LDI. However, the expression of dexterity in this article is a piecewise function: kinematic dexterity is a two-segment function, and dynamic dexterity is a six-segment function. Therefore, it is difficult to apply the method of calculating GDI expression as the objective function in this article.

In this article, two eigenvalue expressions are obtained when calculating the kinematics dexterity, and the kinematics local comprehensive dexterity index (KLCDI) is defined (Eq. (49)). When the KLCDI is closer to 0, the kinematics dexterity of the mechanism is better.

(49) \begin{equation} 0\leq \textit{KLCDI}=\left| \sigma _{\boldsymbol{\max }}-\sigma _{\boldsymbol{\min }}\right| =\left| \frac{k\cos \theta _{1}}{1+\left(k\sin \theta _{1}\right)^{2}}-\frac{k\cos \theta _{2}}{1+\left(k\sin \theta _{2}\right)^{2}}\right| \lt \infty \end{equation}

The KGCDI is further defined as Eq. (50), where w is the area (volume) of the dexterous space. Solving the optimal solution of kinematics dexterity performance can be equivalent to the minimum value of KGCDI.

(50) \begin{equation} 0\leq \textit{KGCDI}=\frac{\sum _{\theta _{1}=-26.6^{\circ}}^{26.6^{\circ}}\sum _{\theta _{2}=-26.6^{\circ}}^{26.6^{\circ}}\textit{KLCDI}}{\int _{w}dw}\lt \infty \end{equation}

6.1.2. Proposal of dynamics global comprehensive dexterity index

In this article, three eigenvalues are solved when calculating the dynamics dexterity. According to the relationship between the eigenvalues, there are six LDI expressions. The dynamics local comprehensive dexterity index (DLCDI) is defined. When the DLCDI is closer to 3, the GIE is closer to the positive circle, and the dynamics dexterity of the mechanism is better.

(51) \begin{equation} 1\lt \textit{DLCDI}=\left| \sigma _{1}+\sigma _{2}+\sigma _{3}\right| =\left| \frac{1+\left(k\sin \theta _{1}\right)^{2}}{k\cos \theta _{1}}+\frac{1+\left(k\sin \theta _{2}\right)^{2}}{k\cos \theta _{2}}+1\right| \lt \infty \end{equation}

The DGCDI is further defined as Eq. (52). Solving the optimal solution of dynamics dexterity performance can be equivalent to solving the DGCDI minimum.

(52) \begin{equation} 1\lt \textit{DGCDI}=\frac{\sum _{\theta _{1}=-26.6^{\circ}}^{26.6^{\circ}}\sum _{\theta _{2}=26.6^{\circ}}^{26.6^{\circ}}\textit{DLCDI}}{\int _{w}dw}\lt \infty \end{equation}

6.1.3. Proposal of GCSI

The evaluation of global stiffness necessitates a multidimensional approach. We contend that a high global average value coupled with a low global range is indicative of superior global stiffness performance. Hence, a GCSI, as depicted in Eq. (53), is proposed.

(53) \begin{equation} GCSI=\frac{\overline{K}-\left(K_{\boldsymbol{\max }}-K_{\boldsymbol{\min }}\right)}{10^{6}} \end{equation}

A larger value of GCSI indicates a better comprehensive stiffness performance of the mechanism.

6.1.4. Establishment of optimization model

The optimization objective function, optimization variables, and variable optimization boundary conditions can be obtained as shown in Eq. (54).

(54) \begin{equation} \left\{\begin{array}{l} min\left\{\begin{array}{l} f_{1}=\textit{KGCDI}\left[k\left(r\mathit{,}\,d_{3}\right)\right]\\ f_{2}=\textit{DGCDI}\left[k\left(r\mathit{,}\,d_{3}\right)\right] \end{array}\right.\\ max\left\{f_{3}=GCSI\left[k\left(r\mathit{,}\,d_{3}\right)\right]\right.\\ \\ s\mathit{.}t\mathit{.}\\ 100\leq r\leq 500\\ 50\leq d_{3}\leq 500 \end{array}\right. \end{equation}

Multiobjective optimization process is developed as shown in Figure 10.

Figure 10. Flowchart of multiobjective optimization.

6.2. Optimization analysis

Considering that the objective function in this article is highly nonlinear and discontinuous in the domain, the NSGA3 is used. Compared with NSGA2, its ability is improved as follows:

1. 1. Improve the ability to deal with large-scale problems.

2. 2. Improve the performance of solving nondominated sorting problems.

3. 3. Improve the maintenance of population diversity to obtain a better solution set.

4. 4. Improve the ability to deal with unevenly distributed objective functions.

The algorithm parameter settings are shown in Table 1. Considering the actual convergence effect and the number of objective functions, the population size is selected to 80 to meet the calculation requirements.

Table 1. Algorithm parameters.

The computer configurations and software versions used for multiobjective optimization are shown in Table 2. The optimization process in this study took less than 5 minutes after many tests. Therefore, to achieve the same or even better computational efficiency, the computer configuration adopted in this study should be satisfied at least.

Table 2. Experimental computer parameter configurations and software versions.

Figure 11. Pareto solution plane and objective function optimization process of NSGA3.

Figure 11(a) is the Pareto optimal solution plane of the NSGA3. With the increase of the number of iterations, the optimal solution gradually converges to one point. Figure 11 (b), (c), and (d) is the optimization process of f ₁, f ₂ and f ₃, respectively. The optimization result of f ₁ is 0.17, the optimization result of f ₂ is 2.68, and the optimization result of f ₃ is 1.10.

Figure 12. The optimized global distribution of kinematics dexterity, dynamics dexterity, and stiffness.

6.3. Optimization results analysis

According to the optimization results of the three objective functions, the values of the optimization variables can be further traced back, which are r = 301.6 and d ₃ = 210.0, respectively. The optimized structural parameters are substituted into the kinematics dexterity, dynamics dexterity, and stiffness calculation model, and the results are shown in Figure 12 (a), (b), and (c). The global kinematics dexterity is above 0.8, and the dexterity of most regions is greater than 0.86. At the same time, the global dynamics dexterity is above 0.66, and the dexterity in most regions is greater than 0.7. the global stiffness is above 1.3 × 10⁶, and the max stiffness is greater than 2.0 × 10⁶.

Meanwhile, KGCDI, DGCDI, and GCSI corresponding to the optimized k value is compared with other k values, and the comparison results are shown in Figure 13.

Figure 13. Comparison of the optimized k value with KGCDI, DGCDI, and GCSI with different k values.

7. Conclusion

In this research, a comprehensive performance analysis and subsequent multiobjective structural optimization of the 4RRR parallel mechanism are undertaken. The core achievements of this study are outlined as follows:

1. New mechanism configuration design: Based on the pure rotation characteristics of 3RRR mechanism, an RRR redundant chain is added, and four chains are arranged at 90° spacing to realize the 3 degrees of freedom decoupling pure rotation of the end of the mechanism. At present, this mechanism is suitable for carrying rollers to roll press complex surfaces. In the future, the mechanism will be used as a decoupled pure rotation module, carrying different tools at the moving platform position to complete different machining tasks.

2. Motion Equivalence and Singularity Analysis: It is established that the 4RRR parallel mechanism, for the purpose of motion analysis, exhibits equivalence to the 3RRR parallel mechanism. Furthermore, the velocity Jacobian matrix is derived under various conditions of releasing distinct degrees of freedom. A singularity analysis is conducted, conclusively demonstrating the absence of singular workspaces within the operational range of the mechanism.

3. Performance Indices Analysis and Multiobjective Optimization: The kinematic dexterity, dynamic dexterity, and stiffness of the mechanism are analyzed. By integrating multiple sets of structural parameters into the analysis framework, the variation patterns of these three critical performance indices are elucidated. To evaluate dexterity performance, two novel indexes—KLCDI and DLCDI—are defined. Additionally, for the purpose of multiobjective optimization, global comprehensive indexes encompassing kinematics (kinematics global comprehensive dexterity index), dynamics (dynamics global comprehensive dexterity index), and stiffness (GCSI) are introduced. Employing the NSGA3 algorithm, a multiobjective optimization process is conducted, ultimately yielding optimal structural parameters (r = 301.6 and d ₃ = 210.0) that enable the mechanism to achieve superior overall performance.