A general method for the manipulability analysis of serial robot manipulators

Abstract

Manipulability quantification and evaluation of robot manipulators plays an important role in robotic mechanism design, motion planning, and control. In this article, a general methodology of using manipulability polytopes to characterize the manipulator’s capability to generate velocities, forces, and accelerations at the end-effector is proposed. This methodology can be used to determine the kinematic, mechanical, and dynamic manipulability of both redundant and nonredundant serial manipulators, where influences of gravity, joint velocities and accelerations, and end-effector payload on the dynamic manipulability are incorporated into the derivation of the manipulability polytopes. Since the manipulability polytope is the feasible region defined by all the joint constraints, the proposed method is suitable for solving the problems including both symmetric and asymmetric joint constraints. Furthermore, a systematic architecture and detailed procedures for performance computation with respect to different types of manipulability of serial robot manipulators are presented. Two typical application examples, respectively, for two and three degree-of-freedom planar manipulators are given to illustrate the correctness and capability of the proposed methodology.

Introduction

Serial robot manipulators have been widely used in manufacturing, aerospace, and logistics industries due to their simple mechanisms, large workspaces, and flexible operations. Performance analysis of manipulators is the prerequisite for the mechanism design, actuator selection and matching, operating style assessment, and motion planning and control. Manipulability measures based on performance analysis provide useful information for the quantification of robot’s inherent properties and for the assessment of robot motion, and hence they are used to construct the objective functions in the mathematical models of robot design and motion planning. Some crucial performance metrics to which should be pay special attention in the process of robot development and use mainly include the kinematic, mechanical, and dynamic capabilities of the robot. The kinematic capability usually refers to the dexterousness which is used to identify the manipulator’s capability to generate velocities in various directions of the task space under the constraints coming from the joint velocity limits. Kinematic capability reflects the velocity transmission ratio of robot manipulator from joint space to task space. The mechanical capability of a robot manipulator for a given configuration is defined as the maximum forces that can be applied by the robot end-effector in any directions of the task space with the restricts of the actuator limits and joint velocity and acceleration limits. The dynamic capability of a robot is the ability of its end-effector to perform accelerations in arbitrary directions of task space, which is used to evaluate the mapping from joint driving space to task space with the consideration of the joint torque and velocity limits and the effects of the gravity and end-effector payload. The aim of this article is to provide a general method applicable to the computation of these three widely used performance metrics.

In the context of manipulability analysis, several approaches have been proposed, such as the traditional or improved manipulability ellipsoids and polytopes. Yoshikawa ¹ proposed the concept “manipulability measure” applicable to both redundant and nonredundant manipulators and established a control algorithm for singularity avoidance of redundant manipulator and obstacle avoidance. Yoshikawa ² also investigated some properties of the manipulability measure and introduced the first definition of the manipulability ellipsoid. Moreover, the total manipulability measure of robotic manipulators was separated into translational and rotational manipulability in Yoshikawa’s research. ³ Vahrenkamp and Asfour ⁴ proposed an extension to the Yoshikawa’s ² manipulability ellipsoid measure. Their extended manipulability measure incorporated the joint limits and the self-distance between the manipulator and other parts of the robot and was used to represent the capability distribution in the robot’s workspace. The aforementioned studies mainly focused on the kinematic manipulability ellipsoid obtained by the application of the singular value decomposition of the Jacobian matrix relating the joint velocities to the end-effector velocities. In fact, the kinematic manipulability ellipsoid approach maps a unit sphere constrained by the normalized joint velocity limits in actuator space to an ellipsoid bounding the end-effector velocity in any direction of task space. Since the manipulability ellipsoid method is inexact in transforming the joint limits into task space, Lee ⁵ used polytope to quantify the motion capability of a robot and proposed a practical method to construct the exact manipulability polytope. Manipulability polytope method has been used for velocity workspace analysis for multiple arm robotic system in Lee’s research, ⁶ in which the problem of the difference and asymmetry of joint limits has been considered and handled. In recent research, Jaquier et al. ⁷ achieved geometry-aware manipulability learning, tracking, and transfer. They established a novel manipulability transfer framework and provided a way by which robots can learn and reproduce manipulability ellipsoid from expert demonstrations. Choi et al. ⁸ proposed a simple desired manipulability ellipsoid used for generating the motion of a manipulator with its trajectory and determining its optimal posture.

Some approaches have also been proposed in the literature to measure the force or mechanical capability of robot manipulators. In Chiu’s work, ⁹ force manipulability ellipsoid was used for measuring the capability of manipulator postures and matching the manipulator’s posture to the task requirements. Chiacchio et al. ¹⁰ proposed an algorithm to correctly evaluate the force manipulability polytope in task space and a new definition of the force manipulability ellipsoid. Their work focused on the static force capability for redundant manipulators. Okabe ¹¹ investigated the translating manipulability polytope of kinematic redundant manipulator and discovered that the translation vector is a function of velocity and acceleration of the extended task space. Nokleby et al. ¹² provided a methodology that allows for the actuator limits and includes both analytical and optimization-based methods to determine the force capability of nonredundantly and redundantly actuated parallel manipulators. Finotello et al. ¹³ used their proposed polytope method in which the translational and rotational analyses are separated to analyze the kinetostatic performances of robot manipulators. Mahjoub and Fahim ¹⁴ took the effect of the gravity into the static behavior analysis of robot manipulators and then integrated the generalized weight that represents the gravity-induced forces acting on individual links of robot with the force ellipsoid to obtain the true force that a manipulator can apply in its task space. Chiacchio et al. ¹⁵ demonstrated that the gravitational forces embedded into the derivation of the ellipsoid can produce a translation of the ellipsoid along all task space directions rather than cause a compression in the ellipsoid volume. Wei and Gao ^16,17 proposed an approach to calculate the output force capacity polytope of multi-degree-of-freedom (DOF) manipulators based on the linear programming method and the geometric method and used their introduced approach to generate the output force capacity polytope and select the driving parameters for a three-DOF excavating manipulator.

When performing kinematic and static manipulability analysis, we neglect the manipulator dynamics. However, in some cases where the robot manipulator must work in the heavy-load, high-speed, and high-precision situations, the manipulator dynamics cannot be completely ignored. Yoshikawa ¹⁸ firstly proposed the concept of dynamic manipulability measure with the consideration of the manipulator dynamics to measure the manipulating ability of the robot end-effector. In Chiacchio and Concilio’s work, ^19,20 a new definition of manipulability ellipsoid that can lead to more correct results and can be used for the case of singular configuration was given to evaluate manipulating capability in terms of end-effector accelerations for redundant manipulators. Rosenstein and Grupen ²¹ demonstrated that the internal motion of a manipulator and its end-effector movement have a complex and nonnegligible effect on dynamic manipulability. Yokokohji et al. ²² extended the dynamic manipulability ellipsoid to the dynamic performance analysis of multifingered grasping systems. Their research proved that the internal forces affect both the volume and the offset of the manipulability ellipsoid bounded the end-effector accelerations, while the gravity forces only lead to an offset of the manipulability ellipsoid. Similar conclusions about the influence of gravity on the dynamic manipulability ellipsoid have also been given by Chiacchio et al. ¹⁵ Zou et al. ²³ used convex polytope method to quantify the dynamic manipulating capability of a hydraulic manipulator in terms of end-effector forces. Zhou et al. ²⁴ conducted the dynamic manipulability analysis using ellipsoid method for multiarm space robot and further extended their method to the closed-loop system composed of the space robot and the captured target. In their research, dynamic manipulability measure was obtained based on dynamic manipulability factor and dynamic manipulability ellipsoid, and the effects of some dynamics parameters, joint variables, and structural geometry parameters on the dynamic manipulability measure were investigated. Chen et al. ²⁵ proposed a dynamic manipulability measure method to assess the unfixed manipulation of an object using multiarm space robot through frictional contacts. Azad et al. ²⁶ studied the effects of the weighting matrix on dynamic manipulability of robots.

In the existing work mentioned above, the three categories of manipulability are analyzed, respectively, by their corresponding methods. There is still no one method that is generally suitable for kinematic, mechanical, and dynamic manipulability. In addition, some dynamics constraints, nonlinear items, and the applicability to both redundant and nonredundant manipulators have not been systematically and adequately taken into account. In this article, we introduce a general method for performance analysis of serial robot manipulators using the manipulability polytope which is a convex feasible region defined by the joint constraints and manipulator kinematics and dynamics. Compared with previous studies, our method is universal for quantifying the end-effector manipulability in terms of velocities, forces, and accelerations in various directions of task space for both redundant and nonredundant serial manipulators. Kinematic manipulability is analyzed by constructing task-space velocity polytope based on kinematic equations, while mechanical and dynamic manipulability analyses are conducted by defining task-space force and acceleration polytopes based on dynamic equations, respectively. When investigating the end-effector accelerations, the influences of Coriolis and centrifugal forces, gravity, and payload are taken into account. In the case of end-effector force capability evaluation, the effects of joint accelerations, gravity, and Coriolis and centrifugal forces are considered. Nonlinear items are handled by linear approximation method, to make it convenient for defining task-space polytopes using feasible region construction method, which is a more generic way for the manipulator with asymmetric joint constraints and redundant DOFs. Moreover, a systematic architecture and detailed procedures are necessary for performance computation with respect to different types of manipulability of serial robot manipulators. Finally, our proposed method is applied to two numerical examples, namely, two- and three-DOF planar manipulators, and it is proven correct and effective.

The remainder of this article includes (1) the basic statement and mathematical models of several manipulator manipulability computation problems; (2) the proposed polytope method for manipulability analysis in terms of task-space velocities, accelerations, and forces; (3) the performance analysis architecture and procedures applicable to any given configuration, end-effector operating trajectories, and entire workspace of a manipulator; (4) two case studies for assessing the performance of two- and three-DOF planar manipulators. Finally, the conclusions are drawn.

Problem statement

The goal of this research is to quantify and characterize the manipulability of serial robot manipulators in task space, including the capabilities of manipulator to generate end-effector velocities, forces, and accelerations in various directions. End-effector output performance is determined mainly by the input joint velocities or driving torques of the joint space, but the end-effector payload, gravity, and some nonlinear terms in dynamic equations also have nonnegligible effects on the mechanical and dynamic capability. Such influence factors arising from driving capability of joint actuators and dynamic behavior of manipulator have to be considered thoroughly and handled effectively to obtain more accurate and realistic manipulator performance.

Robot manipulator transforms the driving velocities or forces provided by the actuators into the task velocities, forces, and accelerations. Thus, as the base of the performance computation, manipulator’s kinematic and dynamic models must be established. For a serial manipulator with n-DOF, its kinematic and dynamic equations can be expressed as

$\dot{\mathbf{p}}=\mathbf{J}(\mathbf{q})\dot{\mathbf{q}}$ 1 $\ddot{\mathbf{p}}=\dot{\mathbf{J}}(\mathbf{q},\dot{\mathbf{q}})\dot{\mathbf{q}}+\mathbf{J}(\mathbf{q})\ddot{\mathbf{q}}$ 2 $M(\mathbf{q})\ddot{\mathbf{q}}+\mathbf{C}(\mathbf{q},\dot{\mathbf{q}})\dot{\mathbf{q}}+J_{(q)}^T{F}_e+G(q)=\tau$ 3

In above equations (1) to (3), $\dot{\mathbf{p}}\in {\Re }^{m\times 1}$ and $\ddot{\mathbf{p}}\in {\Re }^{m\times 1}$ are the vectors of the end-effector velocities and accelerations; $q$ , $\dot{\mathbf{q}}$ , and $\ddot{\mathbf{q}}$ are the n-dimensional vectors of the joint positions, velocities, and accelerations, respectively; $J(q)\in {\Re }^{m\times n}$ is the Jacobian matrix; $J(q)$ can be further differentiated to obtain $\dot{\mathbf{J}}(\mathbf{q},\dot{\mathbf{q}})\in {\Re }^{m\times n}$ ; $M(q)\in {\Re }^{n\times n}$ is the mass matrix; $C(\mathbf{q},\dot{\mathbf{q}})\in {\Re }^{n\times 1}$ represents the Coriolis and centrifugal torques; $F_e\in {\Re }^{m\times 1}$ , $G(q)\in {\Re }^{n\times 1}$ , and $\tau \in {\Re }^{n\times 1}$ are the vectors of the end-effector payload, manipulator gravity, and joint torques, respectively; and m represents the number of the dimensions of the task space and satisfies $m\le n$ .

Combining equations (2) and (3), the end-effector acceleration vector can be further solved as

$\ddot{\mathbf{p}}\text{=}J{M}^{-1}\tau +\dot{\mathbf{J}}\dot{\mathbf{q}}-J{M}^{-1}\mathbf{C}\dot{\mathbf{q}}-J{M}^{-1}{J}^T{F}_e-J{M}^{-1}G$ 4

The forces manipulator that can apply in task space can also be given as

$F_e\text{=}{(J^T)}^{-1}\tau -{(J^T)}^{-1}\mathbf{M}\ddot{\mathbf{q}}-{(J^T)}^{-1}\mathbf{C}\dot{\mathbf{q}}-{(J^T)}^{-1}G$ 5

Considering the constraints of the joint velocity limits, a mathematical computation model of the end-effector velocity in any direction of the task space can be formulated as

$\begin{array}{l}\underset{\dot{\mathbf{q}}}{\max }\hspace{0.17em}{\Vert {\dot{\mathbf{p}}}_{{\varphi }_k}\Vert }_2\hspace{0.17em}\text{for}\hspace{0.17em}{\varphi }_k={\varphi }_0+\frac{2k\pi }{n_{\varphi }}\text{and}\hspace{0.17em}k=0,1,2,\mathrm{...},n_{\varphi }\hspace{0.17em}\\ \text{subject}\hspace{0.17em}\text{to}\\ {\dot{\mathbf{q}}}^L\le \dot{\mathbf{q}}\le {\dot{\mathbf{q}}}^U\end{array}$ 6

where ${\Vert \cdot \Vert }_2$ denotes the Euclidean norm, ${\varphi }_k$ is the direction angle of the end-effector velocity ${\dot{\mathbf{p}}}_{{\varphi }_k}$ , ${\varphi }_0=0$ , and $n_{\varphi }$ is the number of the direction angle steps, the subscripts L and U represent the lower and upper limits of the corresponding vectors.

When investigating the end-effector acceleration performance, an optimization problem subject to maximizing the magnitude of the acceleration can be given as

$\begin{array}{l}\underset{\dot{\mathbf{q}},\tau }{\max }\hspace{0.17em}{\Vert {\ddot{\mathbf{p}}}_{{\varphi }_k}\Vert }_2\hspace{0.17em}\text{for}\hspace{0.17em}{\varphi }_k={\varphi }_0+\frac{2k\pi }{n_{\varphi }}\text{and}\hspace{0.17em}k=0,1,2,\mathrm{...},n_{\varphi }\hspace{0.17em}\\ \text{subject}\hspace{0.17em}\text{to}\\ \{\begin{array}{c}{\dot{\mathbf{q}}}^L\le \dot{\mathbf{q}}\le {\dot{\mathbf{q}}}^U\\ {\tau }^L\le \tau \le {\tau }^U\end{array}\end{array}$ 7

where ${\ddot{\mathbf{p}}}_{{\varphi }_k}$ represents the end-effector acceleration vector with respect to the given direction angle ${\varphi }_k$ .

Similarly, to identify the manipulator’s capability to apply force in any direction of task space, an optimization-based mathematical model can also be formulated as

$\begin{array}{l}\underset{\dot{\mathbf{q}},\ddot{\mathbf{q}},\tau }{\max }\hspace{0.17em}{\Vert F_e{ }_{{\varphi }_k}\Vert }_2\hspace{0.17em}\text{for}\hspace{0.17em}{\varphi }_k={\varphi }_0+\frac{2k\pi }{n_{\varphi }}\text{and}\hspace{0.17em}k=0,1,2,\mathrm{...},n_{\varphi }\hspace{0.17em}\\ \text{subject}\hspace{0.17em}\text{to}\\ \{\begin{array}{c}{\dot{\mathbf{q}}}^L\le \dot{\mathbf{q}}\le {\dot{\mathbf{q}}}^U\\ {\ddot{\mathbf{q}}}^L\le \ddot{\mathbf{q}}\le {\ddot{\mathbf{q}}}^U\\ {\tau }^L\le \tau \le {\tau }^U\end{array}\end{array}$ 8

where $F_e{ }_{{\varphi }_k}$ is the end-effector force vector in the given direction specified by ${\varphi }_k$ .

In the above optimization models, the transformation relationship between the joint space and the end-effector task space is established. It is clear that equations (1), (4) and (5) are the implicit equality constraints of equations (6), (7) and (8), respectively. The optimization-based mathematical models can easily reveal that manipulability analysis for robot manipulators is to calculate the end-effector output capability for velocities, accelerations, and forces in various directions under the joint driving capability limitations and the robot kinematics and dynamics constraints. These models have an advantage on the revelation of the physical significance, but their solving processes are the time-consuming and tedious undertakings. Especially, some challenges caused by the nonlinear terms included in the dynamics constraints cannot be effectively tackled with the existing ellipsoid and polytope methods.

Proposed method

In this section, a systemic and general method to identify the mechanical and motion capability of a serial manipulator in its output space is introduced in detail based on the feasible region construction method of the linear programming problem. A linearization method of the nonlinear terms in the dynamic equations of a manipulator is presented to simplify the manipulability polytope construction. Several typical problems in manipulability polytope construction such as redundant actuation, singular configuration, and asymmetric joint constraints are also discussed.

Feasible region construction method for determining manipulability

To give a better introduction of the application of the feasible region construction method to the manipulability analysis, a methodology of using it to analyze the kinematic and static performance is presented firstly. Kinematic performance analysis is used to determine the end-effector velocity capability under the known joint velocity constraints, as formulated in equations (1) and (6). The feasible region of the end-effector velocity is a vector set of the maximum allowable velocities in various directions of the output space, and it can be defined as

${\Omega }_{\dot{\mathbf{p}}}=\left\{\dot{\mathbf{p}}|{\dot{\mathbf{q}}}^L\le {\mathbf{J}}^{-1}\dot{\mathbf{p}}\le {\dot{\mathbf{q}}}^U\right\}$ 9

In the static performance analysis problem, manipulator capability to apply forces in different directions of end-effector output space is to be found with the constraints of the maximum allowable actuated torques or forces. Compared with the dynamic performance analysis, inertial force and Coriolis and centrifugal forces included in equation (5) are ignored in static force capability calculation. Thus, the feasible region of the end-effector forces manipulator can apply can be described by

${\Omega }_{Fe}\text{=}\left\{F_e|{\tau }^L\le J^T{F}_e\le {\tau }^U\right\}$ 10

In equation (10), the term represents the effect of gravity is not included, and it will be considered in the manner of feasible region translation when conducting the force capability analysis.

According to equation (4), a relation between the end-effector accelerations and the actuated joint torques can be established as

${\Omega }_{\ddot{\mathbf{p}}}\text{=}\left\{\ddot{\mathbf{p}}|{\tau }^L\le M{J}^{-1}\ddot{\mathbf{p}}\le {\tau }^U\right\}$ 11

In equation (11), the terms arise due to the end-effector payloads, Coriolis and gravity, and other inertia forces are not considered to ensure the relation has a linear form. When carrying out the dynamic performance analysis, all those ignored terms will be added.

For a planar robot, its end-effector has a two-dimensional output space, that is, $\dot{\mathbf{p}}\in {\Re }^{2\times 1}$ , $F_e\in {\Re }^{2\times 1}$ , and $\ddot{\mathbf{p}}\in {\Re }^{2\times 1}$ . From equations (9) to (11), it can easily be seen that the feasible regions of the end-effector velocities, forces, and accelerations are bounded by 2n linear inequality constraints coming from the joint velocity and actuated torque limits, respectively. Figure 1 shows the feasible region of a three-DOF planar manipulator. Obviously, the feasible region called the manipulability polytope is the intersection of all feasible regions, respectively, bounded by six linear inequality constraints given in equation (9). In fact, for a three-DOF planar manipulator, the number of the intersection points of any two linear equality constraints selecting from the joint velocity limits is 12. By substituting all these intersection points into the inequalities included in equation (9) and selecting the points satisfying all the inequality constraints, the vertices of the manipulability polytope can be determined. Furthermore, for an n-DOF manipulator whose end-effector acts in an m-dimensional task space, its intersection points of any m linear equality constraints selecting from the joint velocity limits can be computed by solving $C_n^m\cdot 2^m$ linear systems of equations, where $C_n^m$ represents the number of m-combinations from a given set has n elements. Similarly, the obtained $C_n^m\cdot 2^m$ intersection points must also be substituted into the n linear inequality constraints given in equation (9) to identify and select the satisfied vertices. The manipulability polytope in ${\Re }^m$ is bounded by not less than m pairs of hyperplanes, in which each hyperplane is determined by the vertices satisfying all the joint velocity limits.

Figure 1.

The feasible region of the end-effector velocity $\dot{\mathbf{P}}$ .

The polytope constructing method given above can provide an effective mean for manipulator performance visualization, but it cannot directly quantify the manipulator capability in various directions except for those at the vertices of the polytope. Since the manipulator capability to generate velocity in any direction can be characterized by the distance from the origin of the reference frame to the point on the boundary of the polytope in corresponding direction, velocity performance of entire output space can also be obtained by computing the distances in various direction angles ${\varphi }_k$ , $k=0,1,2,\mathrm{...},n_{\varphi }$ , as shown in Figure 2. The range of the direction angle ${\varphi }_k$ is from 0 to $2\pi$ , and the angle step size is $2\pi /n_{\varphi }$ . Corresponding to each direction angle, end-effector velocity satisfies

${\dot{\mathbf{p}}}_{Bi}\hspace{0.17em}=({\dot{p}}_{Bix},\hspace{0.17em}{\dot{p}}_{Biy}),\hspace{0.17em}\frac{{\dot{p}}_{Biy}}{{\dot{p}}_{Bix}}=\tan {\varphi }_k,\hspace{0.17em}i=1,2,\mathrm{...},n_c$ 12

where n_c is the number of the effective constraints in the direction determined by ${\varphi }_k$ . The maximum end-effector velocities bounded by the joint velocity limits in the direction angle ${\varphi }_k$ are represented as ${\dot{\mathbf{p}}}_{Bi}$ , $i=1,2,\mathrm{...},n_c$ , and can be found from

${\overline{\mathbf{J}}}_i{\dot{\mathbf{p}}}_{Bi}\hspace{0.17em}={\dot{q}}_{Bi}$ 13

where ${\dot{q}}_{Bi}$ is the ith joint velocity limit which can be substituted by the upper or lower joint velocity limits, ${\overline{\mathbf{J}}}_i$ is the ith row of $J^{-1}$ . The maximum end-effector velocity on the polytope boundary corresponding to the direction ${\varphi }_k$ is the minimum of ${\dot{\mathbf{p}}}_{Bi}$ , $i=1,2,\mathrm{...},n_c$ ; thus, the polytope boundary is the set of the maximum end-effector velocities in various directions, written as

${\Omega }_{{\dot{\mathbf{p}}}_B}\text{=}\left\{{\dot{\mathbf{p}}}_B{|\Vert {\dot{\mathbf{p}}}_B\Vert }_2=\min \left({\Vert {\dot{\mathbf{p}}}_{B1}\Vert }_2,\mathrm{...},{\Vert {\dot{\mathbf{p}}}_{Bni}\Vert }_2,\mathrm{...},{\Vert {\dot{\mathbf{p}}}_{B{n}_c}\Vert }_2\right),\hspace{0.17em}\frac{{\dot{p}}_{Biy}}{{\dot{p}}_{Bix}}=\tan {\varphi }_k,\hspace{0.17em}{\varphi }_k\in \left[0,2\pi \right]\hspace{0.17em}\right\}$ 14

Figure 2.

Computation of performance in any direction.

Based on equations (10) and (11), the end-effector force and acceleration polytopes only determined by the actuated joint torques can also be constructed by using the above method. However, it is clear from equations (4) and (5) that only the mapping between the end-effector force/acceleration and the actuated joint torques cannot exactly reflect the manipulator output capability. According to equation (5), the end-effector force can be rewritten as

$F_e\text{=}{(J^T)}^{-1}\tau +F_{e-\ddot{q}}+F_{e-C}+F_{e-G}$ 15

Where $F_{e-\ddot{q}}=-{(J^T)}^{-1}\mathbf{M}\ddot{\mathbf{q}}$ , $F_{e-C}=-{(J^T)}^{-1}\mathbf{C}\dot{\mathbf{q}}$ , and $F_{e-G}=-{(J^T)}^{-1}G$ . The terms $F_{e-\ddot{q}}$ , $F_{e-C}$ , and $F_{e-G}$ , respectively, represent the effects of the joint accelerations, the Coriolis and centrifugal torques (or the joint velocities), and the gravity on the end-effector force, and they will change the magnitude of the output forces in various directions of the task space.

For a given manipulator configuration, the effect due to the gravity is a fixed value which will translate the center of the polytope away from the origin where $F_e=0$ , but that produced due to the joint accelerations is different in each direction. Thus, the shape and size of the new polytope with the effect of the joint accelerations will also be changed, compared with that of the original polytope generated only by the actuated joint torque constraints. Considering that $F_{e-\ddot{q}}$ is a linear vector function of the bounded joint acceleration $\ddot{\mathbf{q}}$ , its feasible region is given as

${\Omega }_{F_{e-\ddot{q}}}\text{=}\left\{F_{e-\ddot{q}}|{\ddot{\mathbf{q}}}^L\le -M^{-1}{J}^T{F}_{e-\ddot{q}}\le {\ddot{\mathbf{q}}}^U\right\}$ 16

The set ${\Omega }_{F_{e-\ddot{q}}}$ bounds the offset forces due to the joint acceleration in various directions. To determine the boundary value of ${\Omega }_{F_{e-\ddot{q}}}$ in any direction, the similar method to calculate the boundary of the end-effector velocity polytope as provided in equations (12) to (14) can be used.

Unfortunately, the elements of $F_{e-C}$ are the nonlinear and dependent functions of the joint velocity $\dot{\mathbf{q}}$ ; thus, its value range cannot be computed by directly substituting the joint velocity limits ${\dot{\mathbf{q}}}^L$ and ${\dot{\mathbf{q}}}^U$ for $\dot{\mathbf{q}}$ and constructing the feasible region of the inequalities. Let $\widehat{\mathbf{C}}\text{=}\mathbf{C}\dot{\mathbf{q}}$ and $\widehat{\mathbf{C}}{=[\widehat{\text{c}}}_1\hspace{0.17em}{\widehat{\text{c}}}_2\mathrm{...}{ \widehat{\text{c}}}_n{\text{]}}^T$ , where ${\widehat{\text{c}}}_i$ (i=1,2,…, n) are the elements of the vector $\widehat{\mathbf{C}}$ . It can be seen that each element of $\widehat{\mathbf{C}}$ is a quadratic function in $\dot{\mathbf{q}}$ , as follows

${\widehat{c}}_i={\displaystyle \sum _{k=1}^n{\displaystyle \sum _{l=1}^n\left({\kappa }_{i,kl}\cdot {\dot{q}}_k{\dot{q}}_l\right)}},(i=1,\mathrm{...},n)$ 17

where ${\kappa }_{i,kl}$ is a coefficient of the corresponding quadratic term. In reality, all the elements of $\widehat{\mathbf{C}}$ are nonlinear and dependent on each other. When the norm of $\widehat{\mathbf{C}}$ has the maximum value, the end-effector force polytope determined by the actuated joint torques, the joint accelerations, and the gravity will generate the maximum translating displacement. This translation generated in a specific direction cannot reflect the effect of the joint velocities on the end-effector force capability in various directions. However, to determine the translation displacement in various directions, a number of nonlinear constrained programming problems will lead to challenging computation and lots of time-consuming processes. Naturally, a method transforming the nonlinear and dependent problem into linear and independent problem has to be developed. Detailed theory and process of this method will be introduced in the next subsection. By using this method, the lower and upper bounds of $\widehat{\mathbf{C}}$ can be, respectively, determined as ${\widehat{\mathbf{C}}}^L=[{\widehat{\text{c}}}_1^L\hspace{0.17em}{\widehat{\text{c}}}_2^L\mathrm{...}{\widehat{\text{c}}}_n^L{\text{]}}^T$ and ${\widehat{\mathbf{C}}}^U=[{\widehat{\text{c}}}_1^U\hspace{0.17em}{\widehat{\text{c}}}_2^U\mathrm{...}{\widehat{\text{c}}}_n^U{\text{]}}^T$ . In addition, it is assumed that all elements of ${\widehat{\mathbf{C}}}^L$ and ${\widehat{\mathbf{C}}}^U$ are independent of each other. Obviously, the feasible region of $F_{e-C}$ can be defined as

${\Omega }_{F_{e-C}}\text{=}\left\{F_{e-C}|{\widehat{\mathbf{C}}}^L\le -J^T{F}_{e-C}\le {\widehat{\mathbf{C}}}^U\right\}$ 18

Therefore, the end-effector polytope characterizes that the manipulator’s capability to apply the forces in different directions is the sum of the bias due to the gravity and the polytopes determined, respectively, by the actuated joint torques, joint accelerations, and joint velocities. The values of the points on the boundaries of the above constructed polytopes can be calculated by the similar method used to obtain the boundary values of the end-effector velocity polytope as given in equations (12) to (14).

Here, the acceleration polytope in the output space is also constructed using the similar method. It is assumed that the end-effector payload is determinate when computing the end-effector accelerations. To determine the acceleration feasible region in task space corresponding to the joint actuated torque limits and the joint velocity limits, equation (4) is rearranged as

$\ddot{\mathbf{p}}\text{=}J{M}^{-1}\tau +{\ddot{\mathbf{p}}}_{\dot{q}}+{\ddot{\mathbf{p}}}_{Fe}+{\ddot{\mathbf{p}}}_G$ 19

where ${\ddot{\mathbf{p}}}_{\dot{q}}=\dot{\mathbf{J}}\dot{\mathbf{q}}-J{M}^{-1}\mathbf{C}\dot{\mathbf{q}}$ , ${\ddot{\mathbf{p}}}_{Fe}=-J{M}^{-1}{J}^T{F}_e$ , and ${\ddot{\mathbf{p}}}_G=-J{M}^{-1}G$ . ${\ddot{\mathbf{p}}}_{Fe}$ and ${\ddot{\mathbf{p}}}_G$ represent the translation quantities due, respectively, to the given end-effector payload and the gravity. ${\ddot{\mathbf{p}}}_{\dot{q}}$ represents the translation quantity due to the joint velocities, and its values are different in various directions. The feasible region of ${\ddot{\mathbf{p}}}_{\dot{q}}$ is calculated just like the feasible region of $F_{e-C}$ . The end-effector acceleration ${\ddot{\mathbf{p}}}_{\dot{q}}$ due to the joint velocities can be reformulated as

${\ddot{\mathbf{p}}}_{\dot{q}}=J{M}^{-1}\tilde{\mathbf{C}}$ 20

where $\tilde{\mathbf{C}}=\mathbf{M}{\mathbf{J}}^{-1}\dot{\mathbf{J}}\dot{\mathbf{q}}-\mathbf{C}\dot{\mathbf{q}}$ and $\tilde{\mathbf{C}}$ is a quadratic function in $\dot{\mathbf{q}}$ . Let ${\tilde{\mathbf{C}}}^L=[{\tilde{\text{c}}}_1^L\hspace{0.17em}{\tilde{\text{c}}}_2^L\mathrm{...}{\tilde{\text{c}}}_n^L{\text{]}}^T$ and ${\tilde{\mathbf{C}}}^U=[{\tilde{\text{c}}}_1^U\hspace{0.17em}{\tilde{\text{c}}}_2^U\mathrm{...}{\tilde{\text{c}}}_n^U{\text{]}}^T$ , respectively, be the lower and upper bonds of $\tilde{\mathbf{C}}$ , and the values of their elements will be calculated in detail in the next subsection. Furthermore, the feasible region of ${\ddot{\mathbf{p}}}_{\dot{q}}$ can be expressed as

${\mathbf{\Omega }}_{{\ddot{\mathbf{p}}}_{\dot{q}}}\text{=}\left\{{\ddot{\mathbf{p}}}_{\dot{q}}|{\tilde{\mathbf{C}}}^L\le \mathbf{M}{\mathbf{J}}^{-1}{\ddot{\mathbf{p}}}_{\dot{q}}\le {\tilde{\mathbf{C}}}^U\right\}$ 21

Moreover, the final end-effector acceleration polytope is obtained by translating an overall displacement determined by the end-effector payload and the gravity on the basis of the sum polytope constructed by the actuated joint torque constraints and the joint velocity limits as given in equations (11) and (21).

Linearization method for the joint velocity items in dynamic equations

It has been presented that each element of the vectors representing the effect of the joint velocities is nonlinear and dependent with each other when calculating the end-effector forces and accelerations. For the expressions $F_{e-C}=-{(J^T)}^{-1}\mathbf{C}\dot{\mathbf{q}}$ and ${\ddot{\mathbf{p}}}_{\dot{q}}=\dot{\mathbf{J}}\dot{\mathbf{q}}-J{M}^{-1}\mathbf{C}\dot{\mathbf{q}}$ , considering that both $C$ and $M^{-1}$ are bounded, the boundedness of $F_{e-C}$ and ${\ddot{\mathbf{p}}}_{\dot{q}}$ can be proved. Hence, both $F_{e-C}$ and ${\ddot{\mathbf{p}}}_{\dot{q}}$ are bounded by corresponding quadratic terms in $\dot{\mathbf{q}}$ , as follows ²¹

$\{\begin{array}{c}{\Vert F_{e-C}\Vert }_2\le \alpha (q){\Vert \dot{\mathbf{q}}\Vert }_2^2\\ {\Vert {\ddot{\mathbf{p}}}_{\dot{q}}\Vert }_2\le \beta (q){\Vert \dot{\mathbf{q}}\Vert }_2^2\end{array}$ 22

where both $\alpha (q)$ and $\beta (q)$ are the scalar functions of $q$ . By the upper bounds of $F_{e-C}$ and ${\ddot{\mathbf{p}}}_{\dot{q}}$ specified in equation (22), the determinate offsets corresponding to the end-effector forces and accelerations can be obtained. However, this Euclidean norm-based method cannot reflect the effect of the joint velocities on the end-effector force and acceleration performance in different directions of output space. The magnitude and direction of the translation are fixed in the output space for a given manipulator configuration with the determinate joint velocity limits. To solve this problem, a method to linearize the joint velocity items in dynamic equations and identify its effect on the manipulator end-effector performance in different direction should be proposed.

As equations (18) and (21) expressed, once the lower and upper bounds of $\widehat{\mathbf{C}}$ and $\tilde{\mathbf{C}}$ have been determined, the force and acceleration polytopes identified by the joint velocities can be obtained. Although the results computed by this method are weakly amplified due to the bounds extension of $\widehat{\mathbf{C}}$ and $\tilde{\mathbf{C}}$ in linearization and decoupling process, the obtained polytopes can effectively quantify the effect of joint velocities on the end-effector capability to apply forces and accelerations in various directions and can be easily added to the polytopes generated by the joint actuated torques. Lower and upper bounds of each element of $\widehat{\mathbf{C}}$ and $\tilde{\mathbf{C}}$ must be searched among the stationary points and boundary points of the quadratics of the corresponding elements. The lower and upper bounds of each element of $\widehat{\mathbf{C}}$ and $\tilde{\mathbf{C}}$ are their minimum and maximum values in the allowable ranges of the joint velocities, respectively. Thus, when the bounds of each element of $\widehat{\mathbf{C}}$ and $\tilde{\mathbf{C}}$ are determined, the corresponding extreme points satisfy the following three conditions:

(1) The extreme points are the intersection points of the joint velocity limits. In this case, the extreme points are selected from the combinations of the joint velocity limits. For an n-DOF manipulator, each extreme point is an n-dimensional vector, and each of its elements can be assigned two values: the lower or upper bounds of the corresponding joint velocity. Therefore, the total number of n-combinations from a given set including 2n joint velocity limits is $2^n$ . In all these combinations, the extreme points are those at which the minimum and maximum values of the corresponding element of $\widehat{\mathbf{C}}$ and $\tilde{\mathbf{C}}$ are obtained, respectively.

(2) The extreme points are exactly selected from the stationary points. In this case, for each element of $\widehat{\mathbf{C}}$ and $\tilde{\mathbf{C}}$ , at the stationary points, all the first partial derivatives of the elements ${\widehat{\text{c}}}_i$ and ${\tilde{\text{c}}}_i$ with respect to each joint velocity ${\dot{q}}_j$ are equal to zero, that is, $\partial {\widehat{\text{c}}}_i/\partial {\dot{q}}_j=0$ , $\partial {\tilde{\text{c}}}_i/\partial {\dot{q}}_j=0$ , $(i,j=1,\mathrm{...},n)$ . Among the obtained stationary points, those allowing ${\widehat{\text{c}}}_i$ or ${\tilde{\text{c}}}_i$ to take the minimum or maximum values are the extreme points.

(3) The extreme points are located on the convex or concave boundaries of the allowable region defined by the joint velocity limits. In this case, some elements of the extreme points are the upper or lower bounds of the joint velocities, and the others are those whose first partial derivatives are zero. Assuming that the number of the equations $\partial {\widehat{\text{c}}}_i/\partial {\dot{q}}_j=0$ or $\partial {\tilde{\text{c}}}_i/\partial {\dot{q}}_j=0$ is n_p for any given ${\widehat{\text{c}}}_i$ and ${\tilde{\text{c}}}_i$ and the solutions of these equations belong to the interval constrained by the corresponding joint velocity limits, the number of the elements located on the lower or upper bounds of the joint velocities for the alternative extreme points is $n-n_p$ , then the number of the alternative extreme points is at most $C_n^{n_p}{2}^{n-n_p}$ . Since n_p can take any integer from 0 to n, for any elements ${\widehat{\text{c}}}_i$ and ${\tilde{\text{c}}}_i$ , respectively, of $\widehat{\mathbf{C}}$ and $\tilde{\mathbf{C}}$ , the number of the possible alternative extreme points is at most

$C_n^0{2}^n+C_n^1{2}^{n-1}+C_n^2{2}^{n-2}+\mathrm{...}+C_n^{n_p}{2}^{n-n_p}+\mathrm{...}+C_n^{n-1}{2}^1+C_n^n{2}^0=3^n$ 23

Obviously, this case includes the two cases mentioned previously: in case (1), the number of extreme points is $C_n^0{2}^n$ , and in case (2), $C_n^n{2}^0$ .

Considering $\widehat{\mathbf{C}}\in {\Re }^{n\times 1}$ and $\tilde{\mathbf{C}}\in {\Re }^{n\times 1}$ , the maximum number of the candidate extreme points is $n\hspace{0.17em}\cdot 3^n$ . By comparing the function values of the candidate extreme points, the bounds of each element of $\widehat{\mathbf{C}}$ and $\tilde{\mathbf{C}}$ can be determined. Furthermore, the feasible regions of the end-effector force and acceleration determined by the joint velocities can be constructed according to equations (18) and (21). For the bound vectors ${\widehat{\mathbf{C}}}^L$ , ${\widehat{\mathbf{C}}}^U$ , ${\tilde{\mathbf{C}}}^L$ , and ${\tilde{\mathbf{C}}}^U$ , any two elements of each vector are independent and irrelevant; thus, the processes of the feasible region generation are linear and time-saving.

Analysis of several typical cases

To enhance the generality of the proposed method, several typical problems may encounter in performance analysis, such as the redundancy of the serial manipulator, singular configuration, and asymmetric joint actuated velocity or torque constraints, which are discussed in detail in this section.

Analysis of a case with redundancy

For a redundantly actuated serial manipulator, let n and m, respectively, denote the number of the actuated joints and the task space dimensions, and $n>m$ . In this case, the manipulator Jacobian $J$ is not a square matrix; thus, the inverse Jacobian $J^{-1}$ is replaced by the Moore–Penrose pseudo-inverse of the Jacobian $J$ , denoted $J^{†}$ . It is clear that the terms $M{J}^{†}$ in equations (11) and (16) and $M^{-1}{J}^T$ in equation (16) are $n\times m$ matrices in the polytope construction method outlined above. As presented in previous section, the vertices of the end-effector performance polytopes are obtained by solving $C_n^m\cdot 2^m$ linear systems of equations and selecting those satisfying all the corresponding constraints. Moreover, any point on the boundary of the polytope is determined by the minimum of all the distances truncated by the linear constraints in the corresponding direction denoted by ${\varphi }_k$ , as shown in equation (14) and Figure 2. Whether the joint actuated space is redundant or not, manipulability polytopes for the different end-effector performance can be effectively generated using the proposed feasible region construction method. Therefore, this method can be generalize to the manipulability characterization and evaluation of the serial manipulators, but the result correctness cannot be ensure for the redundant manipulator due to the approximation property of the pseudo-inverse of the Jacobian.

There is also another method to solve the problem of manipulability analysis of redundant manipulator, namely, linear combination or projection-based method, as shown in studies by Lee ⁵ and Wei and Gao. ^16,17 In Lee’s research, manipulability polytope was constructed by linear combination of column vectors of Jacobian matrix or by projecting the joint constraint into a plane made of row vectors of Jacobian matrix. To avoid the incorrectness due to the pseudo-inverse of the Jacobian of redundant manipulator, Lee’s method can be used to construct the end-effector velocity polytope defined by equations (9) and the end-effector acceleration polytopes defined by equations (11) and (21), and the item $J{M}^{-1}\tau$ in equation (19), respectively. For the end-effector force polytopes defined by equations (10), (16), and (18), respectively, our proposed feasible region construction method can give exact results without the computation of pseudo-inverse of the Jacobian. Therefore, our proposed combination method, in which the feasible region construction method is used to generate the end-effector force polytope while the Lee’s method is used to obtain the end-effector velocity and acceleration polytopes, can not only enhance the correctness of the results but also ensure the universality of our method in performance analysis of redundant manipulator.

Analysis of a case with singularity

When a nonredundant manipulator is in a singular configuration, the Jacobian $J$ is a $n\times n$ square matrix, but it is not full rank and $\det \left(J\right)=0$ . When a redundant manipulator is in a singular configuration, $\text{rank}\left(J\right)\le m$ . Since the Jacobian $J$ is not full rank, the matrices $M{J}^{†}$ and $M^{-1}{J}^T$ are also not full rank. Thus, the Moore–Penrose pseudo-inverse formulations of these matrixes should be used to reflect the relationships between the end-effector performance and the joint actuated velocities or torques. Although the manipulator is in the singular configuration, the joint velocity and torque constraints should be satisfied. By the pseudoinverse of the corresponding matrix, the contribution and constraint of each actuated joint to the end-effector performance can be visually and geometrically characterized by using the proposed feasible region construction method. Furthermore, the boundaries of the end-effector performance can also be determined despite the manipulator cannot apply velocities or forces or accelerations in $m-rank\left(J\right)$ principle directions of the task space. Similarly, our proposed combination method of the feasible region construction method and the Lee’s method can also be used to generate the end-effector polytopes of the manipulator in singular configurations to ensure the process generality and results correctness.

Analysis of a case with asymmetric joint constraints

In the traditional methods, whatever the ellipsoid or polytope, the joint limits are normalized in the unitary way even though the lower and upper limits of each joint are asymmetric. In this case, the asymmetric joint limits must be transformed to the formulations which are easy to be described in the normalized way. Then, the augmented matrix or the scaled Jacobian obtained by multiplying the Jacobian matrix with the diagonal matrix whose elements are specified by the arguments is used to replace the Jacobian matrix to develop the relation between the joint space and the task space. Through the augmented matrix, the admission joint actuated velocity or torque spaces described by the unit hyperspheres (for ellipsoid method) or hypercubes (for polytope method) are transformed into the hyperellipsoids or polytopes in task space. In our method, the joint limits are directly substituted into the constraint equations to construct the end-effector feasible output space; thus, it is not necessary to handle the asymmetric joint limits using the unitary way. Similarly, even though the lower and upper bounds of $\widehat{\mathbf{C}}$ and $\tilde{\mathbf{C}}$ are asymmetric, there is no need to normalize them. Obviously, compared to the traditional method, our proposed method has better generality.

Performance analysis architecture and procedures

Based on our proposed general performance analysis method, an overall architecture applicable to manipulability assessment for any given posture, practical or simulated end-effector operating trajectories, and entire workspace of a robot manipulator is established, as shown in Figure 3.

Figure 3.

Architecture for performance analysis.

Performance analysis as illustrated in Figure 3 is performed by a computer program, in which algorithm with respect to each module listed in the architecture is developed. The main procedures for performing evaluation are outlined in the following steps.

Initialization: Give the kinematics and dynamics parameters of the robot manipulator, the end-effector trajectory and manipulator postures, the manipulator motion characteristics, and the external loads acting on the manipulator

Kinematics and dynamics modeling: Establish the manipulator’s kinematics and dynamics models based on the parameters given in the previous step. The manipulator end-effector position can be described by applying the Denavit and Hartenberg (D-H) procedure and the matrix transformation method. Furthermore, the manipulator’s motion equations including kinematic and dynamic equations can be obtained by using the algorithm of Newton–Euler’s formulation.

Position and workspace analysis. The position of the points on the manipulator mechanism is related to the joint angles. The relation between the position of any point and the joint angles is defined based on the D-H coordinate system transformation. Workspace boundary is determined by the maximum reachable range of the manipulator end-effector, which is related to the joint angle limits. Workspace envelope is drawn by position and workspace analysis.

Trajectory, joint space and workspace discretization. The given trajectory is firstly discretized into a series of end-effector working points, then manipulator postures corresponding to each end-effector point belongs to the trajectory are solving. There are two methods to obtain the discretized series of end-effector working points in workspace. One is to obtain the end-effector working points in workspace by discretizing joint angles and direct position analysis, the other is to generate working points through uniformly discretizing the domain enveloped by the workspace boundary. In the former case, the obtained end-effector working points are disorderly and randomly distributed in the workspace, but a one-to-one correspondence between the working points and the postures is established. On the contrary, in the latter case, the obtained working points are uniformly and orderly distributed in the workspace, but each working point corresponds to many or even infinite postures obtained by inverse kinematics analysis.

Inverse kinematics solving: Conduct inverse kinematics analysis for all the working points included in the given trajectory or obtained by workspace discretization. A series of manipulator postures and end-effector orientations corresponding to each working point are obtained to prepare for performance analysis.

General manipulability analysis method programming: Develop computer program modules for general manipulability analysis based on our proposed method and established architecture (Figure 3).

Manipulability computation: Computing manipulability measures and construct output polytopes for the postures which can be arbitrarily specified or obtained by discretizing joint space or solving inverse kinematics for discretized working points in workspace in terms of end-effector velocities, forces and accelerations.

Visualization and statistic analysis of manipulator performance: Visualizing the performance of the manipulator postures which belong to the trajectory or correspond to any working point by manipulability polytopes, and then performing statistic analysis for the manipulability measures of each working points and generating performance spectrum plots for the manipulability measures in the entire workspace.

Case studies

To demonstrate and confirm the correctness and effectiveness of our proposed methodology, two practical applications for performance analysis, respectively, of two- and three-DOF planar manipulators are used in this section. The details on the computation of the matrices and vectors included in the manipulator kinematics and dynamics, which support the performance analysis, can be seen in Appendices 1 and 2.

Case 1: A two-DOF planar manipulator

The parameters used for numerical simulation of a two-DOF manipulator are provided in Table 1. Now, enough foundation is established to evaluate the manipulator’s manipulability to perform. Our objective is to show how to use our proposed method and established architecture as stated in previous sections to analyze the manipulator’s manipulability for any given posture, operating or simulated trajectories and the entire workspace in terms of end-effector velocities, forces and accelerations.

Table 1.

The parameters used for numerical simulation of a two-DOF manipulator.

Parameter and symbol	Value and unit	Parameter and symbol	Value and unit
Length of link 1 (a 1)	1 (m)	Lower bound vector of joint torques ( ${\tau }^L$ )	[−1000, −800] (N/m)
Length of link 2 (a 2)	0.6 (m)	Upper bound vector of joint torques ( ${\tau }^U$ )	[1000, 800] (N/m)
Mass of link 1 (m 1)	20 (kg)	Lower bound vector of joint angles ( $q^L$ )	${[-\frac{\pi }{2},-\frac{4\pi }{9}]}^T$ (rad)
Mass of link 2 (m 2)	12 (kg)	Upper bound vector of joint angles ( $q^U$ )	${[\frac{\pi }{2},\frac{4\pi }{9}]}^T$ (rad)
Moment of inertia of link 1 ( $I_{C_1}$ )	$\frac{m_1{a}_1^2}{12}$ (kg·m2)	Lower bound vector of joint velocities ( ${\dot{\mathbf{q}}}^L$ )	${[-\frac{\pi }{6},-\frac{\pi }{6}]}^T$ (rad/s)
Moment of inertia of link 2 ( $I_{C_2}$ )	$\frac{m_2{a}_2^2}{12}$ (kg·m2)	Upper bound vector of joint velocities ( ${\dot{\mathbf{q}}}^U$ )	${[\frac{\pi }{6},\frac{\pi }{6}]}^T$ (rad/s)

DOF: degree of freedom.

Figures 4 and 5 show the manipulator’s capability to generate velocities in various directions along the horizontal and arc trajectories, respectively. Each end-effector working point on the trajectories corresponds to two manipulator postures, except for those not only on the workspace boundary but also corresponding to the singular postures or those in the domain in which one of two postures corresponding to any working point cannot be achieved due to the joint limits. Comparing the polytopes generated by our proposed method with the traditional ellipsoids, it can be seen that the ellipsoids contained as subdomains of the polytopes underestimate the achievable end-effector velocities. Two manipulator postures corresponding to the same working point have different capability to generate velocities in task space. In Figure 4, both polytopes and ellipsoids become narrow gradually in horizontal direction with the unfolding of the manipulator. Until the joint angle between the two links of the manipulator reaches zero, the polytope together with the ellipsoid degenerate into two segments, respectively. In Figure 5(a), the first type of manipulator posture cannot reach the top domain of the workspace due to the joint limits. Similarly, the second type of manipulator posture cannot reach the bottom domain of the workspace as shown in Figure 5(b). In Figure 5, it can also be seen that the manipulator postures achieve the given arc trajectory have the same capability to apply velocities in various directions of the task space, when the joint angles between the two links of these postures are the same.

Figure 4.

Velocity manipulability ellipsoids (blue) and polytopes (red) of different operating postures for a horizontal trajectory. Workspace boundary (gray).

Figure 5.

Velocity manipulability ellipsoids (blue) and polytopes (red) of different operating postures for an arc trajectory. Workspace boundary (gray).

Figure 6 shows the spectrum plots of two performance measures for two type of manipulator postures. One measure s_v represents the area of the velocity polytope of any manipulator posture, and the other measure r_v is defined by the ratio between the minimum value and the maximum value on the velocity polytope boundary. The former can quantify the global capability of any posture to apply velocity in task space, and the latter can measure the directional uniformity of output capability. The spectrum plots of measures s_v and r_v suggest that performance values are distributed in a radial pattern in the workspace for both two type of manipulator postures. These measure values decrease with increasing operating radius from the left boundary of the reachable workspace to the right. The smaller the operating radius is, the better global output capability and directional uniformity can be achieved.

Figure 6.

Spectrum plots of velocity performance measures for two type of manipulator postures.)

In the force output capability analysis, the lower and upper bounds of the joint accelerations are, respectively, assigned as ${\ddot{\mathbf{q}}}^L={[-\frac{\pi }{3},-\frac{\pi }{3}]}^T$ (rad/s²) and ${\ddot{\mathbf{q}}}^U={[\frac{\pi }{3},\frac{\pi }{3}]}^T$ (rad/s²). Figure 7 shows the componential and resulting ellipsoids and polytopes of output force capability obtained in the manipulator posture $q={[\frac{\pi }{9},-\frac{\pi }{3}]}^T$ , under the influences of the gravity, joint velocities, and joint accelerations. It is quite evident that the effect of the gravity leads to the prominent translations of ellipsoid and polytope, and this effect is greater than that caused by the joint velocities and accelerations. For the given manipulator posture and limits of joint velocities and accelerations, the effects of joint velocities and accelerations are small enough to be ignored. However, for the lightweight and high-speed manipulators, the effect of the gravity will relatively decrease, while the effects of the joint velocity and acceleration will increase to the extent that must be taken into account.

Figure 7.

End-effector force manipulability ellipsoids (a) and polytopes (b) for the posture $q={[\frac{\pi }{9},-\frac{\pi }{3}]}^T$ .

The end-effector force manipulability ellipsoids and polytopes plots given by Figure 8 show the output capability difference between the two types of operating styles for a horizontal trajectory. It is clear that the force capability manipulator can provide in horizontal direction increases gradually with the unfolding of the manipulator, while those in vertical direction decrease slightly. In fact, when the two links of the manipulator are collinear, which is a singular posture, the maximum forces manipulator can apply in a horizontal direction to overcome the external loads no longer depend on the limits of the joint actuating torques but on the manipulator structure strength, at the same time the maximum vertical forces the manipulator can apply are defined by the joint torque limits and the length of the links.

Figure 8.

End-effector force manipulability ellipsoids and polytopes of different operating postures for a horizontal trajectory.

The resulting ellipsoids and polytopes to quantify the force output capability of two type of operating styles for an arc trajectory are plotted in Figure 9. It can clearly be seen that two operating styles can provide different force output capability. The postures with the same joint angles between the two links have the same force output capability in task space. Furthermore, for a given working point, the optimal operating postures can be selected from the two operating styles according to the demands of the objects or task.

Figure 9.

End-effector force manipulability ellipsoids and polytopes of different operating postures for an arc trajectory.

The spectrum plots of force performance measure quantifying the directional uniformity of the force capability for two type of manipulator postures are shown in Figure 10. Similar to the spectrum plots of the velocity performance measure r_v , spectrum plots of the force performance measure r _Fe also extend in a radial pattern in the workspace outward to the right boundary, from the left boundary. It indicates that the greater operating radius is, the worse isotropy of the measure r _Fe for any operating posture is obtained in task space.

Figure 10.

Spectrum plots of force performance measures for two types of manipulator postures.

Figure 11 shows the acceleration manipulability ellipsoids and polytopes corresponding to the manipulator posture $q={[\frac{\pi }{9},-\frac{\pi }{3}]}^T$ with the consideration of the gravity and the joint velocity, where the payload is assigned as zero. It can be seen from the plots that the effect of the gravity makes the ellipsoids or polytopes shift, with no rotation in any direction and no compression in the volume. In addition, the effect of the joint velocity is insignificant for the given manipulator posture.

Figure 11.

End-effector acceleration manipulability ellipsoids (a) and polytopes (b) for the posture $q={[\frac{\pi }{9},-\frac{\pi }{3}]}^T$ .

Figures 12 and 13 show how the manipulator postures affect the end-effector acceleration output capability, respectively, for a horizontal and an arc trajectory. Note that both the acceleration ellipsoids and polytopes are narrow areas in task space. This suggests that acceleration output capability in special directions of task space has an obvious advantage over that in other directions. Some significant characteristics similar to that of velocity and force output capability shown in Figures 4, 5, 8 and 9 can also be seen in Figures 12 and 13; for the sake of simplification, we would not restate them here.

Figure 12.

End-effector acceleration manipulability ellipsoids and polytopes of different operating postures for a horizontal trajectory.

Figure 13.

End-effector acceleration manipulability ellipsoids and polytopes of different operating postures for an arc trajectory.

Figure 14 shows the distribution characteristics of the directional uniformity measure of the acceleration output capability. It is clear that the measure values in the top and bottom of the workspace are greater than that in the middle, and the values are decreased gradually from the left boundary of the reachable workspace to the right.

Figure 14.

Spectrum plots of acceleration performance measures for two types of manipulator postures.

Case 2: A three-DOF planar manipulator

To illustrate our proposed manipulability analysis approach and compare the difference of two type of entire workspace performance evaluation methods, we present an example of a three-DOF planar manipulator in this section. Table 2 provides in detail the parameters used in our example.

Table 2.

The parameters used for numerical simulation of a three-DOF manipulator.

Parameter and symbol	Value and unit	Parameter and symbol	Value and unit
Length of link 1 (a 1)	1 (m)	Lower bound vector of joint torques ( ${\tau }^L$ )	[−400, −200, −200] (N/m)
Length of link 2 (a 2)	0.6 (m)	Upper bound vector of joint torques ( ${\tau }^U$ )	[400, 200, 200] (N/m)
Length of link 3 (a 3)	0.4 (m)	Lower bound vector of joint angles ( $q^L$ )	${[\frac{\pi }{6},-\frac{\pi }{3},-\frac{\pi }{3}]}^T$ (rad)
Mass of link 1 (m 1)	20 (kg)	Upper bound vector of joint angles ( $q^U$ )	${[\frac{5\pi }{6},\frac{\pi }{3},\frac{\pi }{3}]}^T$ (rad)
Mass of link 2 (m 2)	12 (kg)	Lower bound vector of joint velocities ( ${\dot{\mathbf{q}}}^L$ )	${[-\frac{\pi }{6},-\frac{\pi }{6},-\frac{\pi }{6}]}^T$ (rad/s)
Mass of link 3 (m 3)	8 (kg)	Upper bound vector of joint velocities ( ${\dot{\mathbf{q}}}^U$ )	${[\frac{\pi }{6},\frac{\pi }{6},\frac{\pi }{6}]}^T$ (rad/s)
Moment of inertia of link i ( $I_{C_i}$ ), (i = 1, 2, 3)	$\frac{m_i{a}_i^2}{12}$ (kg·m2)

DOF: degree of freedom.

Based on the manipulator geometric and kinematic parameters listed in Table 2, the workspace and working points obtained by two type of discretization methods are shown in Figure 15. In Figure 15(a), the working points in the workspace are obtained by discretizing joint space and direct position analysis for each set of joint angles. It is clear that the working points obtained by this method are disorderly and randomly distributed in the workspace, and each manipulator posture (or each set of joint angles) corresponds to a working point. By contrast, in Figure 15(b), the working points obtained by uniformly discretizing workspace are uniformly and orderly distributed in the workspace, and each working point corresponds to many or even infinite postures obtained by solving the inverse kinematics.

Figure 15.

End-effector working points obtained by two types of discretization methods.

Figure 16 shows the velocity manipulability ellipsoids and polytopes of different operating postures of three given working points. The left and right working points only corresponds to one type of operating postures, while the middle working points can be achieved by two types of postures. It can be noted that the difference of the velocity output capability of the same type of postures for each working point is not so obvious. For the middle working point, two types of operating postures have different predominant directions of output velocities.

Figure 16.

Velocity manipulability ellipsoids and polytopes of different operating postures of three given working points.

Comparison of the directional uniformity measure spectrum plots obtained by discretizing joint space and end-effector workspace is given in Figure 17. Figure 17(a) is the filled contour plot created by interpolating the 3D scattered data containing the directional uniformity measure values on the working points shown in Figure 15(a). Figure 17(b) is directly created based on the directional uniformity measure values on the working points shown in Figure 15(b). As mentioned previously, in Figure 15(a), the working points and manipulator postures satisfy a one-to-one correspondence, while in Figure 15(b) a one-to-many correspondence is satisfied. Therefore, Figure 17(a) is the interpolated results of the directional uniformity measures of the sampled manipulator postures, while Figure 17(b) is filled by the average value of the directional uniformity measures of all the manipulator postures corresponding to every working point. Obviously, two results show a radial distribution pattern except that Figure 17(b) has more clearer and neater contour lines than Figure 17(a).

Figure 17.

Comparison of the directional uniformity measures of velocity manipulability polytopes of different workspace discretization methods.

Figure 18 shows the comparison of the force manipulability ellipsoids and polytopes in the operating posture $q={[\frac{\pi }{3},\frac{\pi }{6},-\frac{\pi }{9}]}^T$ . Figure 19 provides the force manipulability ellipsoids and polytopes of different operating postures of three given working points. Figure 20 compares the distribution characteristics of the directional uniformity measures of force manipulability polytopes obtained, respectively, by joint space discretization and end-effector workspace discretization. In the force manipulability analysis, the lower and upper bounds of the joint accelerations are, respectively, assigned as ${\ddot{\mathbf{q}}}^L={[-\frac{\pi }{3},-\frac{\pi }{3}]}^T$ (rad/s²) and ${\ddot{\mathbf{q}}}^U={[\frac{\pi }{3},\frac{\pi }{3}]}^T$ (rad/s²). From Figures 18 to 20, some conclusions similar to those of two-DOF manipulator and the velocity manipulability of three-DOF manipulator can be drawn. It can clearly be seen that partial boundary of the ellipsoids goes outside of the bounds of the polytopes, and the boundary of the polytopes is no longer a regular shape. This is owing to the influences of the gravity, joint velocity, and acceleration are superimposed on the polytope defined by the joint torques, and this superimposed effect makes the polytopes shift at the same time change their boundary shapes due to the difference of these influences in various directions. However, previous ellipsoid method cannot reflect the difference of each of influences in various directions. In Figure 20, the contour plot obtained by end-effector workspace discretization has a more obviously regional and regular property in spatial distribution than that by joint space discretization.

Figure 18.

End-effector force manipulability ellipsoids (a) and polytopes (b) for the posture $\mathbf{q}={[\frac{\pi }{3},\frac{\pi }{6},-\frac{\pi }{9}]}^T$ .

Figure 19.

Force manipulability ellipsoids and polytopes of different operating postures of three given working points.

Figure 20.

Comparison of the directional uniformity measures of force manipulability polytopes of different workspace discretization methods.

Figure 21 shows the end-effector acceleration manipulability ellipsoids and the polytopes in the operating posture $q={[\frac{\pi }{3},\frac{\pi }{6},-\frac{\pi }{9}]}^T$ and with an end-effector payload of $F_e={[20,60]}^T$ (kg). It can be seen that both the gravity and the end-effector payload make the ellipsoid and polytope shift. However, the shift caused by the payload is not so obvious, since the given payload in this research is small relative to the joint driving capability. In addition, the effect of the joint velocity is smaller than other effects. Figure 22 presents the acceleration manipulability ellipsoids and polytopes of different operating postures of three given working points. Figure 23 shows the contour plots of the directional uniformity measure of the acceleration manipulability polytopes obtained by joint space discretization and end-effector workspace discretization. It is worth noting that the measure values obtained by joint space discretization in Figure 23(a) are much greater than that obtained by end-effector workspace discretization in Figure 23(b). This is because that in Figure 23(a) the contour plot includes the directional uniformity measure values of the all discretized joint postures, while in Figure 23(b) the performance value on each working point represents the average value of the directional uniformity measure values of all operating postures obtained by inverse kinematics analysis. Obviously, the results in Figure 23(a) obtained by as many discretized postures as possible are greater than the average results in Figure 23(b).

Figure 21.

End-effector acceleration manipulability ellipsoids (a) and polytopes (b) for the posture $q={[\frac{\pi }{3},\frac{\pi }{6},-\frac{\pi }{9}]}^T$ .

Figure 22.

Acceleration manipulability ellipsoids and polytopes of different operating postures of three given working points.

Figure 23.

Comparison of the directional uniformity measures of acceleration manipulability polytopes of different workspace discretization methods.

Conclusions

This article introduces a general method for the manipulability analysis of serial robot manipulator based on the construction of the feasible region of the manipulator constraints. Manipulator performance is quantified and characterized in terms of end-effector velocity, force, and acceleration manipulability polytopes. On the one hand, some typical difficulties such as nonlinear constraints, asymmetric joint constraints, redundancy of the manipulator, singularity of Jacobian, and singular configuration are considered and handled in our proposed method. On the other hand, the effects of the gravity, joint velocity and acceleration, and the payload on the end-effector manipulability are also taken into account in the evaluating process. Hence, our proposed method is universal for the analysis of the kinematic, mechanical, and dynamic manipulability of both redundant and nonredundant serial manipulators. Moreover, based on our proposed general performance analysis method, an overall architecture and corresponding procedures applicable to manipulability assessment for any given posture, practical or simulated end-effector operating trajectories, and entire workspace of a robot manipulator are established. The correctness and effectiveness of our proposed methodology are demonstrated by two practical applications, respectively, for two- and three-DOF planar manipulators. Further studies may include extending our method to robot mechanism design and motion planning and optimal control for optimal manipulability.