Kinematical Research of Free-Floating Space-Robot System at Position Level Based on Screw Theory

Abstract

Kinematics of a free-floating space-robot system is a fundamental and complex subject. Problems at the position level, however, are not considered sufficiently because of the nonholonomic property of the system. Current methods cannot handle these problems simply and efficiently. A novel and systematical modeling approach is provided; forward and inverse kinematics at the position level are deduced based on the product of exponentials (POE) formula and conservation of linear momentum. The whole deduction process is concise and clear. More importantly, inertial tensor parameters are not introduced. Then, three situations with different known variables are mainly studied. Due to the complexity of inverse kinematical equations, a numerical method is proposed based on Newton’s iteration method. Two calculation examples are given, a dual-arm planar model and a single-arm spatial model; both forward and inverse kinematical solutions are given, while inverse kinematical results are compared with simulation results of Adams. The results indicate that the proposed methods are quite accurate and efficient.

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

Kinematics is a basic subject of robotics. It includes forward and inverse kinematics at the position and the velocity (differential) level. Mature traditional means have been applied to the research of ground-based robots, such as the Denavit–Hartenberg method [1]. Kinematics of space robots (also referred to as space manipulators) is more complicated due to the base-arm dynamic coupling. For example, the base is completely free in reaction to the arms’ movement when the system is under a free-floating situation [2], which is one of the key researches due to its low fuel consumption. And, motion of the base disturbs the pose (positions and orientations) of the end-effector inversely.

New theories were developed for the kinematics of a space-robot system. Vafa and Dubowsky [3] proposed an ingenious modeling method called the Virtual Manipulator (VM) method. It transforms the real model into a series of virtual-link vectors connected to a fixed base called the Virtual Ground (VG), by using conservation of linear momentum. VG is just located at the centroid of the whole system and is motionless in inertial space. Based on the VM method, Xu et al. [4] presented a concept of “system centroid equivalent manipulator”, whose kinematic model is derived. Then, the trajectory of the balance arm for stabilizing the base centroid is planned. Umetani and Yoshida [5] developed another modeling idea, which defines the generalized Jacobian matrix (GJM). GJM is derived by using the momentum conservation and geometrical relations; it contains inertial parameters of the system. This method is largely used to build the dynamic model of space robots [6–8]. Liang et al. [9] presented another concept called the Dynamically Equivalent Manipulator (DEM). DEM mapped a space manipulator into a fixed-base manipulator which preserves both its dynamic and kinematic properties. It can be physically built for experimental study. Moosavian and Papadopoulos [10] developed and compared two approaches for the kinematics modeling of multiplying robots, called the barycentric vector approach and the direct path method, and concluded that the latter requires significantly less computations. However, inverse kinematics is not discussed.

As angular momentum equations cannot be integrated [11], the space-robot system is nonholonomic. Therefore, the inverse kinematics at the position level cannot be solved only by the pose of the end-effector. This defect leads to that current studies mainly focused on differential kinematics. They usually have to do integration, however, when just joint displacements at the position level are needed in some scenarios, such as path planning [12, 13]. The above-mentioned methods could not solve the kinematics at the position level simply and efficiently. The VM and DEM methods focus on transformation of the model. After transformation, the system is just the same as the ground-based one so that traditional means can be applied. Nevertheless, the transformation process significantly increases the workload, which will be heavier if the system is more complex. The GJM method is derived by using body-fixed vectors, which does not make it easy to program. And its equation is in differential form, which means it cannot be applied to solve the kinematics at the position level directly. What is more, the equations of GJM and DEM contain dynamic parameters, which are not necessary to solve kinematical problems at the position level. Liang and Xu [14] discussed the possibility of an inverse kinematical solution at the position level when there are some constraints and known variables. Joint angles of a single-arm space robot are calculated by using the D-H method. A multiarm space-robot system is not concerned. To this end, this paper focused on the kinematics of a free-floating space-robot system at the position level, aimed at developing a simple and efficient method. Moreover, the method should cover forward and inverse kinematics and single-arm and multiarm space-robot systems.

Screw theory is an efficient alternative to analysis kinematics of a rigid body and mechanism. It was first systematized by R.S. Ball, then employed and improved for study on ground-based manipulators [15, 16] and parallel robots [17, 18]. A screw depicts a dual vector, so it can represent both translation and rotation, which gives screw theory the advantages of unity and simplicity. As a basic and conventional tool in screw theory to study robotic kinematics, the product of exponentials (POE) formula is derived by the successive screw displacement of a kinematic chain. Rocha et al. [19] compared the screw-based method with the D-H method in detail and concluded that the former has advantages over the D-H one, although the D-H method is more commonly used. Nevertheless, few studies on a space-robot system using screw theory are found in current research. Liu and Wu [20] deduced the dynamic model of a space-robot system based on the screw theory and Kane equation. Liu and Wu [21] provided a simple and efficient method to deduce the Jacobian matrix of a free-floating branching robot system, by dividing it into two parts, a ground robot with a fixed base and a locking joint robot with a floating base. Kinematics at the position level is not considered.

Taking the advantages of the screw theory, this paper provides a novel and systematical approach. Forward and inverse kinematical equations of a free-floating space-robot system at the position level are established based on the POE formula and conservation of linear momentum, without employing inertial tensors. And several situations with different known variables are discussed: only position of the base is known, only orientation of the base is known, and both the position and orientation of the base are measurable. A numerical means based on Newton’s iteration method is given, considering that the closed-form inverse solution may be hard to calculate directly when the system is complex. The entire approach can be applied to space robots at free-floating mode without regard to the number of arms.

2. The Model and Basic Theories

2.1. The Model of the Space-Robot System

A space-robot system consists of the base (satellite) and arms (also called manipulators, single or multiple) typically. Each arm consists of several links, joints, and an end-effector (EE). The following assumptions are made to build a concise and accurate model.

(1)

The position and orientation of the base is not controlled, and there are no external forces on the system, which means the system is at free-floating mode

(2)

Each arm is a single chain connected to the base, and the parallel manipulator is not considered

(3)

All links and joints are rigid

(4)

All joints are revolute, and all the rotation angles are limited between $-{180}^{\circ }$ and ${180}^{\circ }$

(5)

The number of dofs of each arm is no greater than six. For a planar model, it is no greater than three. A redundant system is not considered

A model is provided in Figure 1. The definitions of each symbol are given as follows.

(i)

${\sum }_{\mathbf{I}}$ represents the inertial frame and is located at the centroid of the whole system

(ii)

${\sum }_{\text{B}}$ represents the base-fixed frame

(iii)

Vector ${{}^{\mathbf{I}}\mathbf{p}}_0$ represents the centroidal position of the base in the inertial frame

(iv)

Vector ${{}^{\mathbf{I}}\mathbf{p}}_i^k$ represents the centroidal position of the $k^{th}$ arm’s $i^{th}$ link in the inertial frame

(v)

Vector ${{}^{\mathbf{I}}\mathbf{p}}_E^k$ represents the position of the $k^{th}$ arm’s EE in the inertial frame

[figure omitted; refer to PDF]

2.2. Product of the Exponential Formula

The POE formula describes the kinematical relationship by using the exponential form of screws. Only two frames are needed, the tool frame $\left\{T\right\}$ and the referential frame $\left\{S\right\}$. The formula is written as $$\begin{array}{}\text{(1)}& {}_T{}^{S}g\left(\mathbf{\theta }\right)=e^{{\theta }_1{\widehat{\mathit{\xi }}}_1}{e}^{{\theta }_2{\widehat{\mathit{\xi }}}_2}\text{L}{e}^{{\theta }_n{\widehat{\mathit{\xi }}}_n}{}_T{}^{S}g\left(\mathbf{0}\right),\end{array}$$where ${}_T{}^{S}g\left(\mathbf{\theta }\right)$ and ${}_T{}^{S}g\left(\mathbf{0}\right)$ represent the final and initial poses of frame $\left\{T\right\}$ relative to frame $\left\{S\right\}$, respectively. $e^{{\theta }_i{\widehat{\mathit{\xi }}}_i}$ is just the exponential form of rigid motion. ${\mathit{\xi }}_i$ defines the twist of the corresponding joint. For a revolute joint, $$\begin{array}{}\text{(2)}& {\mathbf{\xi }}_i={\left({\mathbf{\omega }}_i,{\mathbf{r}}_i\times {\mathbf{\omega }}_i\right)}^T,\end{array}$$where ${\mathbf{\omega }}_i$ is the unit vector along the axis of the twist and ${\mathbf{r}}_i$ is the position vector of an arbitrary point on the axis.

${}_T{}^{S}g$ can also be expressed as $$\begin{array}{}\text{(3)}& {}_T{}^{S}g=\left({}_T{}^S\text{R},{}^S\mathbf{p}\right)\in \text{SE}\left(3\right),\end{array}$$where ${}_T{}^S\text{R}$ is the rotational matrix and ${}^S\mathbf{p}$ is the position vector.

${}_T{}^{S}g$ satisfies an important rule: $$\begin{array}{}\text{(4)}& {}_T{}^{A}g={}_S{}^{A}g{}_T{}^{S}g.\end{array}$$

It means the poses of frame $\left\{T\right\}$ relative to frame $\left\{S\right\}$ and frame $\left\{A\right\}$ can be transformed mutually.

3. Kinematics Deduction Base on Screw Theory

3.1. Forward Kinematics

In this section, forward kinematical equations of the above-mentioned model are established.

According to the POE formula, the final pose of the centroidal frame of the $k^{th}$ arm’s $i^{th}$ link in the base-fixed frame can be given by $$\begin{array}{}\text{(5)}& {}_{C_i^k}{}^{B}g=e^{{\theta }_1^k{\widehat{\mathit{\xi }}}_1^k}{e}^{{\theta }_2^k{\widehat{\mathit{\xi }}}_2^k}\text{L}{e}^{{\theta }_i^k{\widehat{\mathit{\xi }}}_i^k}{}_{C_i^k}{}^{B}g\left(0\right)=\left[\begin{array}{cc}{{{}^{B}R}_i}^k& {{{}^B\mathbf{p}}_i}^k\\ \mathbf{0}& 1\end{array}\right],\hspace{1em}i=1,2,K,n^k,k=1,2,\cdots ,N,\end{array}$$and the final pose of the $k^{th}$ arm’s EE can be given by $$\begin{array}{}\text{(6)}& {}_{E^k}{}^{B}g=e^{{\theta }_1^k{\widehat{\mathit{\xi }}}_1^k}{e}^{{\theta }_2^k{\widehat{\mathit{\xi }}}_2^k}\cdots e^{{\theta }_n^k{\widehat{\mathit{\xi }}}_n^k}{}_{E^k}{}^{B}g\left(0\right)=\left[\begin{array}{cc}{{{}^{B}R}_E}^k& {{{}^B\mathbf{p}}_E}^k\\ \mathbf{0}& 1\end{array}\right].\end{array}$$

The final pose of the base and centroidal frame of the $k^{th}$ arm’s $i^{th}$ link in the inertial frame can be expressed as $$\begin{array}{c}\text{(7)}& {}_B{}^{I}g=\left[\begin{array}{cc}{{}^{I}R}_0& {{}^I\mathbf{p}}_0\\ \mathbf{0}& 1\end{array}\right],{}_{{C_i}^k}{}^{I}g=\left[\begin{array}{cc}{{{}^{I}R}_i}^k& {{{}^I\mathbf{p}}_i}^k\\ \mathbf{0}& 1\end{array}\right],\end{array}$$where ${{}^{I}R}_0$ is determined by the orientation of the base. If we define the attitude angles (roll, pitch, and yaw) of the base as $\mathbf{\Psi }={\left(\delta ,\beta ,\gamma \right)}^T$, then $$\begin{array}{}\text{(8)}& {{}^{I}R}_0=\left[\begin{array}{ccc}\text{c}\beta \text{c}\gamma & \text{s}\delta \text{s}\beta \text{c}\gamma -\text{c}\delta \text{s}\gamma & \text{c}\delta \text{s}\beta \text{c}\gamma +\text{s}\delta \text{s}\gamma \\ \text{c}\beta \text{s}\gamma & \text{s}\delta \text{s}\beta \text{s}\gamma +\text{c}\delta \text{c}\gamma & \text{c}\delta \text{s}\beta \text{s}\gamma -\text{s}\delta \text{c}\gamma \\ -\text{s}\beta & \text{s}\delta \text{c}\beta & \text{c}\delta \text{c}\beta \end{array}\right],\end{array}$$where $\text{c}$ means $\cos$ and $\text{s}$ means $\sin$.

According to the transformation rule of rigid motion, there is $$\begin{array}{}\text{(9)}& {}_{{C_i}^k}{}^{I}g={}_B{}^{I}g{}_{{C_i}^k}{}^{B}g=\left[\begin{array}{cc}{{}^{I}R}_0& {{}^I\mathbf{p}}_0\\ \mathbf{0}& 1\end{array}\right]\left[\begin{array}{cc}{{{}^{B}R}_i}^k& {{{}^B\mathbf{p}}_i}^k\\ \mathbf{0}& 1\end{array}\right]=\left[\begin{array}{cc}{{}^{I}R}_0{{{}^{B}R}_i}^k& {{}^{I}R}_0{{{}^B\mathbf{p}}_i}^k+{{}^I\mathbf{p}}_0\\ \mathbf{0}& 1\end{array}\right].\end{array}$$

So, $$\begin{array}{}\text{(10)}& \left\{\begin{array}{c}{{{}^{I}R}_i}^k={{}^{I}R}_0{{{}^{B}R}_i}^k,\hfill \\ {{}^I\mathbf{p}}_i={{}^{I}R}_0{{{}^B\mathbf{p}}_i}^k+{{}^I\mathbf{p}}_0,\hfill \end{array}\right.\end{array}$$and the final pose of the $k^{th}$ arm’s EE in the inertial frame can be expressed as $$\begin{array}{}\text{(11)}& {}_{E^k}{}^{I}g={}_B{}^{I}g{}_{E^k}{}^{B}g=\left[\begin{array}{cc}{{}^{I}R}_0{{}^{B}R}_E^k& {{}^{I}R}_0{{}^B\mathbf{p}}_E^k+{{}^I\mathbf{p}}_0\\ \mathbf{0}& 1\end{array}\right].\end{array}$$

In the inertial frame, the centroidal formula of the system is written as $$\begin{array}{}\text{(12)}& m_0{{}^I\mathbf{p}}_0+\sum _{k=1}^N\sum _{i=1}^{n^k}{m_i}^k{{{}^I\mathbf{p}}_i}^k=M{{}^I\mathbf{p}}_C=0.\end{array}$$

Substituting equation (10) into equation (12) yields $$\begin{array}{}\text{(13)}& {{}^I\mathbf{p}}_0=-\frac{{{}^{I}R}_0}{M}\sum _{k=1}^N\sum _{i=1}^{n^k}{m_i}^k{{{}^B\mathbf{p}}_i}^k.\end{array}$$

Based on these deductions, ${}_{E^k}{}^{I}g$ can be calculated while the known variables are changing.

Situation 1.

When both the position and orientation of the base are measurable, equation (11) can be used directly.

Situation 2.

When only the orientation of the base can be obtained, substituting equation (13) into equation (11), there is $$\begin{array}{}\text{(14)}& {}_{E^k}{}^{I}g=\left[\begin{array}{cc}{{}^{I}R}_0{{}^{B}R}_E^k& {{}^{I}R}_0{{}^B\mathbf{p}}_E^k-\frac{{{}^{I}R}_0}{M}\sum _{k=1}^N\sum _{i=1}^{n^k}{m_i}^k{{{}^B\mathbf{p}}_i}^k\\ \mathbf{0}& 1\end{array}\right],\end{array}$$where ${{{}^B\mathbf{p}}_i}^k$ can be calculated by equation (5).

Situation 3.

When only the position of the base can be obtained, substituting equation (8) into equation (13), there is $$\begin{array}{}\text{(15)}& {{}^I\mathbf{p}}_0=-\left[\begin{array}{ccc}\text{c}\beta \text{c}\gamma & \text{s}\delta \text{s}\beta \text{c}\gamma -\text{c}\delta \text{s}\gamma & \text{c}\delta \text{s}\beta \text{c}\gamma +\text{s}\delta \text{s}\gamma \\ \text{c}\beta \text{s}\gamma & \text{s}\delta \text{s}\beta \text{s}\gamma +\text{c}\delta \text{c}\gamma & \text{c}\delta \text{s}\beta \text{s}\gamma -\text{s}\delta \text{c}\gamma \\ -\text{s}\beta & \text{s}\delta \text{c}\beta & \text{c}\delta \text{c}\beta \end{array}\right]\frac{1}{M}\sum _{k=1}^N\sum _{i=1}^{n^k}{m_i}^k{{{}^B\mathbf{p}}_i}^k.\end{array}$$

Let ${{}^I\mathbf{p}}_0={\left[p_{01},p_{02},p_{03}\right]}^T,\left(1/M\right){\sum }_{k=1}^N{\sum }_{i=1}^{n^k}{m_i}^k{{{}^B\mathbf{p}}_i}^k={\mathbf{p}}_t={\left[p_{t1},p_{t2},p_{t3}\right]}^T$, then equation (15) can be written as $$\begin{array}{}\text{(16)}& \left\{\begin{array}{c}\text{c}\beta \text{c}\gamma p_{t1}+\left(\text{s}\delta \text{s}\beta \text{c}\gamma -\text{c}\delta \text{s}\gamma \right)p_{t2}+\left(\text{c}\delta \text{s}\beta \text{c}\gamma +\text{s}\delta \text{s}\gamma \right)p_{t3}+p_{01}=0,\hfill \\ \text{c}\beta \text{s}\gamma p_{t1}+\left(\text{s}\delta \text{s}\beta \text{s}\gamma +\text{c}\delta \text{c}\gamma \right)p_{t2}+\left(\text{c}\delta \text{s}\beta \text{s}\gamma -\text{s}\delta \text{c}\gamma \right)p_{t3}+p_{02}=0,\hfill \\ -\text{s}\beta p_{t1}+\text{s}\delta \text{c}\beta p_{t2}+\text{c}\delta \text{c}\beta p_{t3}+p_{03}=0.\hfill \end{array}\right.\end{array}$$

$\mathbf{\Psi }=\left(\delta ,\beta ,\gamma \right)$ can be calculated using this equation, then ${}_{E^k}{}^{I}g$ can be calculated by equation (11).

3.2. Inverse Kinematics

Firstly, the feasibility of the inverse solution should be demonstrated. The number of state variables of the system is $\left(6+{\sum }_{k=1}^N{n}^k\right)$, including the position and orientation of the base. The number of equations is $6N$, which is equal to the total number of state parameters of all the EEs. The formula of conservation of linear momentum is a holonomic constraint, which can eliminate three independent state variables of the system. If the other three state variables are known, the number of unknown state variables will be no greater than the number of equations $\left({\sum }_{k=1}^N{n}^k\le 6N\right)$. Then, the equations will have limited numbers of solutions. For the planar model, the explanation is similar. Accordingly, the following situations with different known variables are discussed.

Situation 1.

When both the position and the orientation of the base are measurable, substituting equation (6) into equation (11), there is $$\begin{array}{}\text{(17)}& {}_{E^k}{}^{I}g={}_B{}^{I}g{}_{E^k}{}^{B}g={}_B{}^{I}g{e}^{{\theta }_1^k{\widehat{\mathit{\xi }}}_1^k}{e}^{{\theta }_2^k{\widehat{\mathit{\xi }}}_2^k}\text{L}{e}^{{\theta }_n^k{\widehat{\mathit{\xi }}}_n^k}{}_{E^k}{}^{B}g\left(0\right).\end{array}$$

Let, $$\begin{array}{}\text{(18)}& G=e^{{\theta }_1^k{\widehat{\mathit{\xi }}}_1^k}{e}^{{\theta }_2^k{\widehat{\mathit{\xi }}}_2^k}\text{L}{e}^{{\theta }_n^k{\widehat{\mathit{\xi }}}_n^k}={{}_B{}^{I}g}^{-1}{}_{E^k}{}^{I}g{}_{E^k}{}^{B}g{\left(0\right)}^{-1}.\end{array}$$

Joint angles of the $k^{th}$ arm ${\mathbf{\theta }}^k=\left\{{\theta }_i^k\right\}$ can be solved by using the Paden-Kahan subproblem method. It will not be discussed further because it belongs to a maturely studied subject.

Situation 2.

When only the orientation of the base can be obtained, equation (14) can be converted into the following form: $$\begin{array}{}\text{(19)}& \left[\begin{array}{cc}{{}^{I}R}_0{{}^{B}R}_E^k& {{}^{I}R}_0{{}^B\mathbf{p}}_E^k-\frac{{{}^{I}R}_0}{M}\sum _{k=1}^N\sum _{i=1}^{n^k}{m_i}^k{{{}^B\mathbf{p}}_i}^k\\ \mathbf{0}& 1\end{array}\right]-{}_{E^k}{}^{I}g=Y^k\left(\mathbf{\theta }\right)-{}_{E^k}{}^{I}g=0,\end{array}$$where $Y^k\in \text{SE}\left(3\right)$. Twelve elements of $Y^k$ are dependent on the joint angles. Considering that nine elements of the rotational matrix are dependent mutually, let $$\begin{array}{c}\text{(20)}& Y^k=\left[\begin{array}{cccc}{\text{y}}_{11}^k& {\text{y}}_{12}^k& {\text{y}}_{13}^k& {\text{y}}_{14}^k\\ {\text{y}}_{21}^k& {\text{y}}_{22}^k& {\text{y}}_{23}^k& {\text{y}}_{24}^k\\ {\text{y}}_{31}^k& {\text{y}}_{32}^k& {\text{y}}_{33}^k& {\text{y}}_{34}^k\\ & \mathbf{0}& & 1\end{array}\right],{}_{E^k}{}^{I}g=\left[\begin{array}{cccc}{t}_{11}^k& t_{12}^k& t_{13}^k& t_{14}^k\\ t_{21}^k& t_{22}^k& t_{23}^k& t_{24}^k\\ t_{31}^k& t_{32}^k& t_{33}^k& t_{34}^k\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

By using the unit quaternion method, equation (19) can be simplified to $$\begin{array}{c}\text{(21)}& \left\{\begin{array}{c}{f}_1^k\left(\mathbf{\theta }\right)=\frac{1}{2}\sqrt{\left(1+{\text{y}}_{11}^k+{\text{y}}_{22}^k+{\text{y}}_{33}^k\right)}-\frac{1}{2}\sqrt{\left(1+t_{11}^k+t_{22}^k+t_{33}^k\right)}=0,\hfill \\ f_2^k\left(\mathbf{\theta }\right)=\left({\text{y}}_{32}^k-{\text{y}}_{23}^k\right)-\left(t_{32}^k-t_{23}^k\right)=0,\hfill \\ f_3^k\left(\mathbf{\theta }\right)=\left({\text{y}}_{13}^k-{\text{y}}_{31}^k\right)-\left(t_{13}^k-t_{31}^k\right)=0,\hfill \\ f_4^k\left(\mathbf{\theta }\right)=\left({\text{y}}_{21}^k-{\text{y}}_{12}^k\right)-\left(t_{21}^k-t_{12}^k\right)=0,\hfill \end{array}\right.\left\{\begin{array}{c}{f}_5^k\left(\mathbf{\theta }\right)={\text{y}}_{14}^k-t_{14}^k=0,\hfill \\ f_6^k\left(\mathbf{\theta }\right)={\text{y}}_{24}^k-t_{24}^k=0,\hfill \\ f_7^k\left(\mathbf{\theta }\right)={\text{y}}_{34}^k-t_{34}^k=0.\hfill \end{array}\right.\end{array}$$

All the equations of all the arms compose the following equation set: $$\begin{array}{}\text{(22)}& F\left(\mathbf{\theta }\right)=\mathbf{0};\hspace{1em}\mathbf{\theta }=\left\{{{\theta }_i}^k\right\},i=1,2,K,n^k,k=1,2,\cdots ,N.\end{array}$$

All joint angles $\mathbf{\theta }$ can be obtained by solving this equation set. However, it may be hard to get a closed-form solution due to its complexity. A numerical approach is provided based on Newton’s iteration method [22]. Firstly, the Jacobian matrix of the equation set is given by $$\begin{array}{}\text{(23)}& J_F\left(\mathbf{\theta }\right)=\left[\begin{array}{cccccc}\frac{\partial f_1^1}{\partial {\theta }_1^1}& \frac{\partial f_1^1}{\partial {\theta }_2^1}& \cdots & \frac{\partial f_1^1}{\partial {\theta }_1^2}& \begin{array}{c}\cdots \end{array}& \frac{\partial f_1^1}{\partial {\theta }_{n^N}^N}\\ \frac{\partial f_2^1}{\partial {\theta }_1^1}& \frac{\partial f_2^1}{\partial {\theta }_2^1}& \cdots & \frac{\partial f_2^1}{\partial {\theta }_1^2}& \cdots & \frac{\partial f_2^1}{\partial {\theta }_{n^N}^N}\\ ⋮& ⋮& \ddots & ⋮& \ddots & ⋮\\ \frac{\partial f_1^2}{\partial {\theta }_1^1}& \frac{\partial f_1^2}{\partial {\theta }_2^1}& \cdots & \frac{\partial f_1^2}{\partial {\theta }_1^2}& \cdots & \frac{\partial f_1^2}{\partial {\theta }_{n^N}^N}\\ ⋮& ⋮& \ddots & ⋮& \ddots & ⋮\\ \frac{\partial f_7^N}{\partial {\theta }_1^1}& \frac{\partial f_7^N}{\partial {\theta }_2^1}& \cdots & \frac{\partial f_7^N}{\partial {\theta }_1^2}& \cdots & \frac{\partial f_7^N}{\partial {\theta }_{n^N}^N}\end{array}\right].\end{array}$$

$J_F\left(\mathbf{\theta }\right)\in R^{\left(7N\right)\times N_{\text{total}}}$, where $N_{\text{total}}={\sum }_{k=1}^N{n}^k$. Then, the iterative formula can be written as $$\begin{array}{}\text{(24)}& {\mathbf{\theta }}^{i+1}={\mathbf{\theta }}^i-{J_F}^{+}F\left(\mathbf{\theta }\right),\end{array}$$where ${J_F}^{+}$ is the pseudoinverse matrix of $J_F$.

Situation 3.

When only the position of the base can be obtained, the problem is more complex. Firstly, to convert equation (9) into the following form: $$\begin{array}{}\text{(25)}& \left[\begin{array}{cc}{{}^{I}R}_0{{}^{B}R}_E^k& {{}^{I}R}_0{{}^B\mathbf{p}}_E^k-\frac{{{}^{I}R}_0}{M}\sum _{k=1}^N\sum _{i=1}^{n^k}{m_i}^k{{{}^B\mathbf{p}}_i}^k\\ \mathbf{0}& 1\end{array}\right]-{}_{E^k}{}^{I}g=Z^k\left(\mathbf{\theta },\mathbf{\Psi }\right)-{}_{E^k}{}^{I}g=\mathbf{0},\end{array}$$let, $$\begin{array}{}\text{(26)}& Z^k=\left[\begin{array}{cccc}{z}_{11}^k& z_{12}^k& z_{13}^k& z_{14}^k\\ z_{21}^k& z_{22}^k& z_{23}^k& z_{24}^k\\ z_{31}^k& z_{32}^k& z_{33}^k& z_{34}^k\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

Similarly, equation (25) can be simplified to $$\begin{array}{c}\text{(27)}& \left\{\begin{array}{c}{g}_1^k\left(\mathbf{\theta },\mathbf{\Psi }\right)=\frac{1}{2}\sqrt{\left(1+z_{11}^k+z_{22}^k+z_{33}^k\right)}-\frac{1}{2}\sqrt{\left(1+t_{11}^k+t_{22}^k+t_{33}^k\right)}=0,\hfill \\ g_2^k\left(\mathbf{\theta },\mathbf{\Psi }\right)=\left(z_{32}^k-z_{23}^k\right)-\left(t_{32}^k-t_{23}^k\right)=0,\hfill \\ g_3^k\left(\mathbf{\theta },\mathbf{\Psi }\right)=\left(z_{13}^k-z_{31}^k\right)-\left(t_{13}^k-t_{31}^k\right)=0,\hfill \\ g_4^k\left(\mathbf{\theta },\mathbf{\Psi }\right)=\left(z_{21}^k-z_{12}^k\right)-\left(t_{21}^k-t_{12}^k\right)=0,\hfill \end{array}\right.\left\{\begin{array}{c}{g}_5^k\left(\mathbf{\theta },\mathbf{\Psi }\right)=z_{14}^k-t_{14}^k=0,\hfill \\ g_6^k\left(\mathbf{\theta },\mathbf{\Psi }\right)=z_{24}^k-t_{24}^k=0,\hfill \\ g_7^k\left(\mathbf{\theta },\mathbf{\Psi }\right)=z_{34}^k-t_{34}^k=0.\hfill \end{array}\right.\end{array}$$

Equation (16) can be converted into the following form: $$\begin{array}{}\text{(28)}& \left\{\begin{array}{c}{\text{g}}_1^{\prime }\left(\mathbf{\theta },\mathbf{\Psi }\right)=\text{c}\beta \text{c}\gamma p_{t1}+\left(\text{s}\delta \text{s}\beta \text{c}\gamma -\text{c}\delta \text{s}\gamma \right)p_{t2}+\left(\text{c}\delta \text{s}\beta \text{c}\gamma +\text{s}\delta \text{s}\gamma \right)p_{t3}+p_{01}=0,\hfill \\ {\text{g}}_2^{\prime }\left(\mathbf{\theta },\mathbf{\Psi }\right)=\text{c}\beta \text{s}\gamma p_{t1}+\left(\text{s}\delta \text{s}\beta \text{s}\gamma +\text{c}\delta \text{c}\gamma \right)p_{t2}+\left(\text{c}\delta \text{s}\beta \text{s}\gamma -\text{s}\delta \text{c}\gamma \right)p_{t3}+p_{02}=0,\hfill \\ {\text{g}}_3^{\prime }\left(\mathbf{\theta },\mathbf{\Psi }\right)=-\text{s}\beta p_{t1}+\text{s}\delta \text{c}\beta p_{t2}+\text{c}\delta \text{c}\beta p_{t3}+p_{03}=0.\hfill \end{array}\right.\end{array}$$

All the equations of all the arms and equation (28) compose another equation set: $$\begin{array}{}\text{(29)}& G\left(\mathbf{\phi }\right)=\mathbf{0};\hspace{1em}\mathbf{\phi }=\left\{{{\theta }_i}^k,\delta ,\beta ,\gamma \right\},i=1,2,K,n^k,k=1,2,\cdots ,N.\end{array}$$

The Jacobian matrix of the equation set is given by $$\begin{array}{}\text{(30)}& J_G\left(\mathbf{\phi }\right)=\left[\begin{array}{ccccccccc}\frac{\partial g_1^1}{\partial {\theta }_1^1}& \frac{\partial g_1^1}{\partial {\theta }_2^1}& \cdots & \frac{\partial g_1^1}{\partial {\theta }_1^2}& \cdots & \frac{\partial g_1^1}{\partial {\theta }_{n^N}^N}& \frac{\partial g_1^1}{\partial \delta }& \frac{\partial g_1^1}{\partial \beta }& \frac{\partial g_1^1}{\partial \gamma }\\ \frac{\partial g_2^1}{\partial {\theta }_1^1}& \frac{\partial g_2^1}{\partial {\theta }_2^1}& \cdots & \frac{\partial g_2^1}{\partial {\theta }_1^2}& \cdots & \frac{\partial g_2^1}{\partial {\theta }_{n^N}^N}& \frac{\partial g_2^1}{\partial \delta }& \frac{\partial g_2^1}{\partial \beta }& \frac{\partial g_2^1}{\partial \gamma }\\ ⋮& ⋮& \ddots & ⋮& \ddots & ⋮& ⋮& ⋮& ⋮\\ \frac{\partial g_1^2}{\partial {\theta }_1^1}& \frac{\partial g_1^2}{\partial {\theta }_2^1}& \cdots & \frac{\partial g_1^2}{\partial {\theta }_1^2}& \cdots & \frac{\partial g_1^2}{\partial {\theta }_{n^N}^N}& \frac{\partial g_1^2}{\partial \delta }& \frac{\partial g_1^2}{\partial \beta }& \frac{\partial g_1^2}{\partial \gamma }\\ ⋮& ⋮& \ddots & ⋮& \ddots & ⋮& ⋮& ⋮& ⋮\\ \frac{\partial g_7^N}{\partial {\theta }_1^1}& \frac{\partial g_7^N}{\partial {\theta }_2^1}& \cdots & \frac{\partial g_7^N}{\partial {\theta }_1^2}& \cdots & \frac{\partial g_7^N}{\partial {\theta }_{n^N}^N}& \frac{\partial g_7^N}{\partial \delta }& \frac{\partial g_7^N}{\partial \beta }& \frac{\partial g_7^N}{\partial \gamma }\\ \frac{\partial {\text{g}}_1^{\prime }}{\partial {\theta }_1^1}& \frac{\partial {\text{g}}_1^{\prime }}{\partial {\theta }_2^1}& \cdots & \frac{\partial {\text{g}}_1^{\prime }}{\partial {\theta }_1^2}& \cdots & \frac{\partial {\text{g}}_1^{\prime }}{\partial {\theta }_{n^N}^N}& \frac{\partial {\text{g}}_1^{\prime }}{\partial \delta }& \frac{\partial {\text{g}}_1^{\prime }}{\partial \beta }& \frac{\partial {\text{g}}_1^{\prime }}{\partial \gamma }\\ \frac{\partial {\text{g}}_2^{\prime }}{\partial {\theta }_1^1}& \frac{\partial {\text{g}}_2^{\prime }}{\partial {\theta }_2^1}& \cdots & \frac{\partial {\text{g}}_2^{\prime }}{\partial {\theta }_1^2}& \cdots & \frac{\partial {\text{g}}_2^{\prime }}{\partial {\theta }_{n^N}^N}& \frac{\partial {\text{g}}_2^{\prime }}{\partial \delta }& \frac{\partial {\text{g}}_2^{\prime }}{\partial \beta }& \frac{\partial {\text{g}}_2^{\prime }}{\partial \gamma }\\ \frac{\partial {\text{g}}_3^{\prime }}{\partial {\theta }_1^1}& \frac{\partial {\text{g}}_3^{\prime }}{\partial {\theta }_2^1}& \cdots & \frac{\partial {\text{g}}_3^{\prime }}{\partial {\theta }_1^2}& \cdots & \frac{\partial {\text{g}}_3^{\prime }}{\partial {\theta }_{n^N}^N}& \frac{\partial {\text{g}}_3^{\prime }}{\partial \delta }& \frac{\partial {\text{g}}_3^{\prime }}{\partial \beta }& \frac{\partial {\text{g}}_3^{\prime }}{\partial \gamma }\end{array}\right].\end{array}$$

$J_G\left(\mathbf{\phi }\right)\in R^{\left(7N+3\right)\times \left(N_{\text{total}}+3\right)}$. The iterative formula is written as $$\begin{array}{}\text{(31)}& {\mathbf{\phi }}^{i+1}={\mathbf{\phi }}^i-{J_G}^{+}G\left(\mathbf{\phi }\right),\end{array}$$where ${J_G}^{+}$ is the pseudoinverse matrix of $J_G$.

In particular, all the equations in this section can be simplified into kinematical equations of a single-arm space-robot system when $N=1$.

4. Examples and Discussion

4.1. Solution of the Planar Model

Firstly, a simulation example of a three-link dual-arm planar robot is given in this part. The model at the initial configuration is shown in Figure 2.

[figure omitted; refer to PDF]

To make a comparison between the simulation results of MATLAB and Adams, a simplified model of the robot is built in Adams and a certain motion process (from the initial configuration to the final configuration) is run. The final configuration of the Adams model is shown in Figure 3.

[figure omitted; refer to PDF]

Detailed parameters of the model are given as follows: $$\begin{array}{c}\text{(32)}& {\mathbf{\omega }}_1^1={\mathbf{\omega }}_2^1={\mathbf{\omega }}_3^1={\mathbf{\omega }}_1^2={\mathbf{\omega }}_2^2={\mathbf{\omega }}_3^2={\left(0,0,1\right)}^T,{{\mathbf{r}}_1}^1={\left(1.2,0.7,0\right)}^T,\\ {{\mathbf{r}}_2}^1={\left(3.0,0.7,0\right)}^T,\\ {{\mathbf{r}}_3}^1={\left(4.8,0.7,0\right)}^T\\ {{\mathbf{r}}_1}^2={\left(-1.2,0.7,0\right)}^T,\\ {{\mathbf{r}}_2}^2={\left(-3.0,0.7,0\right)}^T,\\ {{\mathbf{r}}_3}^2={\left(-4.8,0.7,0\right)}^T.\end{array}$$

So, $$\begin{array}{c}\text{(33)}& {{\mathbf{\xi }}_1}^1={\left({{\mathbf{\omega }}_1}^1,{{\mathbf{r}}_1}^1\times {{\mathbf{\omega }}_1}^1\right)}^T={\left(0,0,1,0.7,-1.2,0\right)}^T,{{\mathbf{\xi }}_2}^1={\left({{\mathbf{\omega }}_2}^1,{{\mathbf{r}}_2}^1\times {{\mathbf{\omega }}_2}^1\right)}^T={\left(0,0,1,0.7,-3,0\right)}^T,\\ {{\mathbf{\xi }}_3}^1={\left({{\mathbf{\omega }}_3}^1,{{\mathbf{r}}_3}^1\times {{\mathbf{\omega }}_3}^1\right)}^T={\left(0,0,1,0.7,-4.8,0\right)}^T,\\ {{\mathbf{\xi }}_1}^2={\left({{\mathbf{\omega }}_1}^2,{{\mathbf{r}}_1}^2\times {{\mathbf{\omega }}_1}^2\right)}^T={\left(0,0,1,0.7,1.2,0\right)}^T,\\ {{\mathbf{\xi }}_2}^2={\left({{\mathbf{\omega }}_2}^2,{{\mathbf{r}}_2}^2\times {{\mathbf{\omega }}_2}^2\right)}^T={\left(0,0,1,0.7,3,0\right)}^T,\\ {{\mathbf{\xi }}_3}^2={\left({{\mathbf{\omega }}_3}^2,{{\mathbf{r}}_3}^2\times {{\mathbf{\omega }}_3}^2\right)}^T={\left(0,0,1,0.7,4.8,0\right)}^T,\end{array}$$and $m_0=400,{m_1}^1=18,{m_2}^1=18,{m_3}^1=24,{m_1}^2=18,{m_2}^2=18,{m_3}^2=24$, and $M=520$. $$\begin{array}{c}\text{(34)}& {}_{C_1^1}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& 2.1\\ E_{3\times 3}& 0.7\\ & 0\\ \mathbf{0}& 1\end{array}\right],{}_{C_2^1}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& 3.9\\ E_{3\times 3}& 0.7\\ & 0\\ \mathbf{0}& 1\end{array}\right],\\ {}_{C_3^1}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& 5.7\\ E_{3\times 3}& 0.7\\ & 0\\ \mathbf{0}& 1\end{array}\right],\\ {}_{E^1}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& 6.6\\ E_{3\times 3}& 0.7\\ & 0\\ \mathbf{0}& 1\end{array}\right],\\ {}_{C_1^2}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& -2.1\\ E_{3\times 3}& 0.7\\ & 0\\ \mathbf{0}& 1\end{array}\right],\\ {}_{C_2^2}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& -3.9\\ E_{3\times 3}& 0.7\\ & 0\\ \mathbf{0}& 1\end{array}\right],\\ {}_{C_3^2}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& -5.7\\ E_{3\times 3}& 0.7\\ & 0\\ \mathbf{0}& 1\end{array}\right],\\ {}_{E^2}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& -6.6\\ E_{3\times 3}& 0.7\\ & 0\\ \mathbf{0}& 1\end{array}\right].\end{array}$$

4.1.1. Forward Kinematical Solution of the Planar Model

Base on the method in Section 3.1, forward kinematical solution of this planar model is given.

Final joint angles are $\mathbf{\theta }={\left({83.15}^{°},{18.94}^{°},{128.4}^{°},-{37.19}^{°},-{29.44}^{°},-{51.09}^{°}\right)}^T$. Then, according to equation (5), there are $$\begin{array}{c}\text{(35)}& {}_{C_1^1}{}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_1^1& {{}^B\mathbf{p}}_1^1\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_1^1{\widehat{\mathit{\xi }}}_1^1}{}_{C_1^1}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}0.1193& -0.9929& 0& 1.3073\\ 0.9929& 0.1193& 0& 1.5936\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],{}_{C_2^1}{}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_2^1& {{}^B\mathbf{p}}_2^1\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_2^1{\widehat{\mathit{\xi }}}_2^1}{e}^{{\theta }_2^1{\widehat{\mathit{\xi }}}_2^1}{}_{C_2^1}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}-0.2094& -0.9778& 0& 1.2262\\ 0.9778& -0.2094& 0& 1.3672\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],\\ {}_{C_3^1}{}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_3^1& {{}^B\mathbf{p}}_3^1\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_3^1{\widehat{\mathit{\xi }}}_3^1}{e}^{{\theta }_3^1{\widehat{\mathit{\xi }}}_3^1}{e}^{{\theta }_3^1{\widehat{\mathit{\xi }}}_3^1}{}_{C_{3^1}}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}-0.6362& 0.7715& 0& 0.4647\\ -0.7715& -0.6362& 0& 3.5511\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],\\ {}_{C_1^2}{}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_1^2& {{}^B\mathbf{p}}_1^2\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_1^2{\widehat{\mathit{\xi }}}_1^2}{}_{C_1^2}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}0.7966& 0.6045& 0& -1.9170\\ -0.6045& 0.7966& 0& 1.2440\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],\\ {}_{C_2^2}{}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_2^2& {{}^B\mathbf{p}}_2^2\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_2^2{\widehat{\mathit{\xi }}}_2^2}{e}^{{\theta }_2^2{\widehat{\mathit{\xi }}}_2^2}{}_{C_2^2}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}0.3967& 0.9180& 0& -2.9909\\ -0.9180& 0.3967& 0& 2.6142\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],\\ {}_{C_3^2}{}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_3^2& {{}^B\mathbf{p}}_3^2\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_3^2{\widehat{\mathit{\xi }}}_3^2}{e}^{{\theta }_3^2{\widehat{\mathit{\xi }}}_3^2}{e}^{{\theta }_3^2{\widehat{\mathit{\xi }}}_3^2}{}_{C_3^2}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}-0.4652& 0.8852& 0& -2.9293\\ -0.8852& -0.4652& 0& 4.2371\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

According to equation (6), there are $$\begin{array}{c}\text{(36)}& {}_{E^1}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_E^1& {{}^B\mathbf{p}}_E^1\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_1^1{\widehat{\mathit{\xi }}}_1^1}{e}^{{\theta }_2^1{\widehat{\mathit{\xi }}}_2^1}{e}^{{\theta }_3^1{\widehat{\mathit{\xi }}}_3^1}{}_{E^1}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}-0.6362& 0.7715& 0& -0.1079\\ -0.7715& -0.6362& 0& 2.8568\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],{}_{E^2}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_E^2& {{}^B\mathbf{p}}_E^2\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_1^2{\widehat{\mathit{\xi }}}_1^2}{e}^{{\theta }_2^2{\widehat{\mathit{\xi }}}_2^2}{e}^{{\theta }_3^2{\widehat{\mathit{\xi }}}_3^2}{}_{E^2}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}‐0.4652& 0.8852& 0& ‐2.5107\\ -0.8852& ‐0.4652& 0& 5.0338\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

Situation 1.

Both the position and the orientation of the base are known: $$\begin{array}{c}\text{(37)}& {{}^I\mathbf{p}}_0={\left(0.1470,-0.6772,0\right)}^T,\mathbf{\Psi }={\left(0,0,-{4.18}^{\circ }\right)}^T.\end{array}$$

So, $$\begin{array}{}\text{(38)}& {}_B{}^{I}g=\left[\begin{array}{cccc}0.9973& 0.0729& 0& 0.1470\\ -0.0729& 0.9973& 0& -0.6772\\ 0& 0& 1& 0\\ & 0& & 1\end{array}\right].\end{array}$$

According to equation (11), there are $$\begin{array}{c}\text{(39)}& {}_{E^1}{}^{I}g={}_B{}^{I}g{}_{E^1}{}^{B}g=\left[\begin{array}{cccc}-0.6908& 0.7231& 0& 0.2476\\ -0.7231& -0.6908& 0& 2.1798\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],{}_{E^2}{}^{I}g={}_B{}^{I}g{}_{E^2}{}^{B}g=\left[\begin{array}{cccc}-0.5284& 0.8490& 0& -1.9901\\ -0.8490& -0.5284& 0& 4.5261\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

Situation 2.

Only the orientation of the base is obtained: $\mathbf{\Psi }={\left(0,0,-{4.18}^{\circ }\right)}^T$.

So, $$\begin{array}{}\text{(40)}& {{}^{I}R}_0=\left[\begin{array}{ccc}0.9973& 0.0729& 0\\ -0.0729& 0.9973& 0\\ 0& 0& 1\end{array}\right].\end{array}$$

According to equation (14), there are $$\begin{array}{c}\text{(41)}& {}_{E^1}{}^{I}g=\left[\begin{array}{cc}{{}^{I}R}_0{{}^{B}R}_E^1& {{}^{I}R}_0{{}^B\mathbf{p}}_E^1-\frac{{{}^{I}R}_0}{M}\sum _{k=1}^2\sum _{i=1}^3{m}_i^k{{}^B\mathbf{p}}_i^k\\ \mathbf{0}& 1\end{array}\right]=\left[\begin{array}{cccc}-0.6908& 0.7231& 0& 0.2476\\ -0.7231& -0.6908& 0& 2.1798\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],{}_{E^2}{}^{I}g=\left[\begin{array}{cc}{{}^{I}R}_0{{}^{B}R}_E^2& {{}^{I}R}_0{{}^B\mathbf{p}}_E^2-\frac{{{}^{I}R}_0}{M}\sum _{k=1}^2\sum _{i=1}^3{m}_i^k{{}^B\mathbf{p}}_i^k\\ \mathbf{0}& 1\end{array}\right]=\left[\begin{array}{cccc}-0.5284& 0.8490& 0& -1.9901\\ -0.8490& -0.5284& 0& 4.5261\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

Situation 3.

Only the position of the base is obtained: ${{}^I\mathbf{p}}_0={\left(0.1470,-0.6772,0\right)}^T$.

Based on the known parameters, there is $$\begin{array}{}\text{(42)}& {\mathbf{p}}_t=\frac{1}{M}\sum _{k=1}^2\sum _{i=1}^3{m}_i^k{{}^B\mathbf{p}}_i^k={\left(-0.1959,0.6647,0\right)}^T.\end{array}$$

According to equation (16), there is $$\begin{array}{}\text{(43)}& \left\{\begin{array}{c}-0.1959\text{c}\gamma -0.6647\text{s}\gamma +0.1470=0,\hfill \\ -0.1959\text{s}\gamma +0.6647\text{c}\gamma -0.6772=0.\hfill \end{array}\right.\end{array}$$

So, $\gamma =-0.0729=-{4.17}^{\circ }$. Then, $$\begin{array}{}\text{(44)}& {}_B{}^{I}g=\left[\begin{array}{cccc}0.9974& 0.0727& 0& 0.1470\\ -0.0727& 0.9974& 0& -0.6772\\ 0& 0& 1& 0\\ & 0& & 1\end{array}\right].\end{array}$$

According to equation (11), there are $$\begin{array}{c}\text{(45)}& {}_{E^1}{}^{I}g={}_B{}^{I}g{}_{E^1}{}^{B}g=\left[\begin{array}{cccc}-0.6906& 0.7232& 0& 0.2472\\ -0.7232& -0.6906& 0& 2.1798\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],{}_{E^2}{}^{I}g={}_B{}^{I}g{}_{E^2}{}^{B}g=\left[\begin{array}{cccc}‐0.5283& 0.8491& 0& ‐1.9909\\ ‐0.8491& ‐0.5283& 0& 4.5258\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

Comparing the results of these situations, which are almost equal to each other, it can be concluded that the method is feasible and effective. The result of Situation 3 has a subtle error due to the truncation error of angles during the calculation.

4.1.2. Inverse Kinematical Solution of the Planar Model

Based on the method in Section 3.2, iterative programs for each situation are written with MATLAB and are run separately. Firstly, the final pose of the two arms is given by $$\begin{array}{c}\text{(46)}& {}_{E^1}{}^{I}g=\left[\begin{array}{cccc}-0.6908& 0.7231& 0& 0.2476\\ -0.7231& -0.6908& 0& 2.1798\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],{}_{E^2}{}^{I}g=\left[\begin{array}{cccc}-0.5284& 0.8490& 0& -1.9901\\ -0.8490& -0.5284& 0& 4.5261\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

The process of the iterative programs is as follow.

Step 1.

Set the calculating parameters and known variables. Max number of iterations is $N=100$. Convergence tolerance is $\epsilon =0.001$. Known parameters are $\mathbf{\Psi }={\left(0,0,-{4.18}^{\circ }\right)}^T$ (for Situation 2). ${{}^I\mathbf{p}}_0={\left(0.1469,-0.6773,0\right)}^T$ (for Situation 3).

Step 2.

Set initial joint angles. ${\mathbf{\theta }}_0={\left({60}^{°},{30}^{°},{80}^{°},-{60}^{°},-{30}^{°},-{80}^{°}\right)}^T$. For the program of Situation 3, the initial attitude angles of the base should also be given: ${\mathbf{\Psi }}_0={\left(0,0,-{10}^{°}\right)}^T$.

Step 3.

Calculate the Jacobian matrix of each iterative step using equation (23) (for Situation 2) or equation (30) (for Situation 3).

Step 4.

Calculate the angles of each iterative step using equation (24) (for Situation 2) or equation (31) (for Situation 3).

Step 5.

Calculate the residual error and judge whether the iteration is over. If yes, go to the next step. If no, go to step 2.

Step 6.

Calculate the final pose of each arm using the iterative result of the angles.

The results of the Adams simulation and the iterations are listed in Table 1.

Table 1

Simulation results of the Adams simulation and the iterations for the planar model.

Method	Joint angles (°)	Base’s attitude angles	Max error
Adams	$${\left(83.15,18.94,128.4,-37.19,-29.44,-51.09\right)}^T$$	$${\left(0,0,-{4.18}^{°}\right)}^T$$	—
Iteration for Situation 2	$${\left(83.096,19.014,128.351,-37.185,-29.449,-51.105\right)}^T$$	—	0.39%
Iteration for Situation 3	$${\left(83.089,19.020,128.364,-37.214,-29.443,-51.077\right)}^T$$	$${\left(0,0,-{4.163}^{°}\right)}^T$$	0.42%

The comparison shows that the results of the iterations are quite approximate to the real solution, which demonstrates the accuracy of the proposed methods. Actually, comparison between forward kinematical solutions and results of the iterations could also verify the accuracy.

Furthermore, efficiency of the methods is studied. Iterative time and number of iterative steps under certain calculating conditions are recorded while the initial parameters of the iteration of Situation 2 are changing. The results are presented in Table 2.

Table 2

Simulation results of iterative time and number of iterative steps.

Initial joint angles (°)	Iterative time (s)	Steps
$${\left({100}^{°},{50}^{°},{60}^{°},-{80}^{°},-{60}^{°},-{100}^{°}\right)}^T$$	5.94	11
$${\left({60}^{°},{30}^{°},{80}^{°},-{60}^{°},-{30}^{°},-{80}^{°}\right)}^T$$	3.41	6
$${\left({65}^{°},{25}^{°},{95}^{°},-{55}^{°},-{30}^{°},-{60}^{°}\right)}^T$$	3.51	6
$${\left({70}^{°},{23}^{°},{105}^{°},-{45}^{°},-{30}^{°},-{55}^{°}\right)}^T$$	2.97	5
$${\left({75}^{°},{20}^{°},{115}^{°},-{40}^{°},-{30}^{°},-{50}^{°}\right)}^T$$	2.47	4

The results show that the iterative method is fairly efficient, and it could perform better when the initial parameters are closer to the real values. This is determined by the inherent property of Newton’s iterative method.

4.2. Solution of the Spatial Model

To verify the application of the proposed method to a tridimensional model, a calculating example of a single-arm three-link space robot is presented. The Adams model at the initial configuration and the final configuration is shown in Figures 4 and 5, respectively.

[figure omitted; refer to PDF] [figure omitted; refer to PDF]

The parameters of this spatial model are given as follows: $$\begin{array}{c}\text{(47)}& {\mathbf{\omega }}_1={\mathbf{\omega }}_3={\left(0,0,1\right)}^T,{\mathbf{\omega }}_2={\left(0,-1,0\right)}^T\\ {\mathbf{r}}_1={\left(1.2030,0.6929,0\right)}^T,\\ {\mathbf{r}}_2={\left(2.903,0.6929,0\right)}^T,\\ {\mathbf{r}}_3={\left(4.603,0.6929,0\right)}^T\end{array}$$

So, $$\begin{array}{c}\text{(48)}& {\mathbf{\xi }}_1={\left({\mathbf{\omega }}_1,{\mathbf{r}}_1\times {\mathbf{\omega }}_1\right)}^T={\left(0,0,1,0.6929,-1.2030,0\right)}^T,{\mathbf{\xi }}_2={\left({\mathbf{\omega }}_2,{\mathbf{r}}_2\times {\mathbf{\omega }}_2\right)}^T={\left(0,-1,0,0,0,-2.9030\right)}^T,\\ {\mathbf{\xi }}_3={\left({\mathbf{\omega }}_3,{\mathbf{r}}_3\times {\mathbf{\omega }}_3\right)}^T={\left(0,0,1,0.6929,-4.603,0\right)}^T,\end{array}$$and $m_0=400,m_1=15,m_2=15,m_3=24$, and $M=454$. $$\begin{array}{c}\text{(49)}& {}_{C_1}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& 2.275\\ E_{3\times 3}& 0.6929\\ & 0\\ \mathbf{0}& 1\end{array}\right],{}_{C_2}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& 3.975\\ E_{3\times 3}& 0.6929\\ & 0\\ \mathbf{0}& 1\end{array}\right],\\ {}_{C_3}{}^{B}g\left(0\right)=\left[\begin{array}{cc}& 5.675\\ E_{3\times 3}& 0.6929\\ & 0\\ \mathbf{0}& 1\end{array}\right],\\ {}_E^{B}g\left(0\right)=\left[\begin{array}{cc}& 6.306\\ E_{3\times 3}& 0.6829\\ & 0\\ \mathbf{0}& 1\end{array}\right].\end{array}$$

4.2.1. Forward Kinematical Solution of the Spatial Model

The final joint angles are $\mathbf{\theta }={\left(-{67.51}^{°},{88.35}^{°},{21.62}^{°}\right)}^T$. According to equation (5), there are $$\begin{array}{c}\text{(50)}& {}_{C_1}{}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_1& {{}^B\mathbf{p}}_1\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_1{\widehat{\mathit{\xi }}}_1}{}_{C_1}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}0.3825& 0.9239& 0& 1.6131\\ -0.9239& 0.3825& 0& -0.2976\\ 0& 0& 1& 0\\ & \mathbf{0}& & 1\end{array}\right],{}_{C_2}{}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_2& {{}^B\mathbf{p}}_2\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_1{\widehat{\mathit{\xi }}}_1}{e}^{{\theta }_2{\widehat{\mathit{\xi }}}_2}{}_{C_2}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}0.0110& 0.9239& -0.3824& 1.8651\\ -0.0266& 0.3825& 0.9236& -0.9063\\ 0.9996& 0& 0.0288& 1.0716\\ & \mathbf{0}& & 1\end{array}\right],\\ {}_{C_3}{}^{B}g=\left[\begin{array}{cc}{{}^{B}R}_3& {{}^B\mathbf{p}}_3\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_1{\widehat{\mathit{\xi }}}_1}{e}^{{\theta }_2{\widehat{\mathit{\xi }}}_2}{e}^{{\theta }_3{\widehat{\mathit{\xi }}}_3}{}_{C_3}{}^{B}g\left(0\right)=\left[\begin{array}{cccc}0.3507& 0.8549& -0.3824& 2.2479\\ 0.1162& 0.3654& 0.9236& -0.7985\\ 0.9293& -0.3683& 0.0288& 2.6955\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

According to equation (6), there is $$\begin{array}{}\text{(51)}& {}_E{}^{B}g=\left[\begin{array}{cc}{{{}^{B}R}_E}^k& {{{}^B\mathbf{p}}_E}^k\\ \mathbf{0}& 1\end{array}\right]=e^{{\theta }_1{\widehat{\mathit{\xi }}}_1}{e}^{{\theta }_2{\widehat{\mathit{\xi }}}_2}{e}^{{\theta }_3{\widehat{\mathit{\xi }}}_3}{}_E{}^{B}g\left(0\right)=\left[\begin{array}{cccc}0.3507& 0.8549& -0.3824& 2.4602\\ 0.1162& 0.3654& 0.9236& -0.7290\\ 0.9293& -0.3683& 0.0288& 3.2857\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

Situation 1.

Both the position and the orientation of the base are known: ${{}^I\mathbf{p}}_0={\left(-0.3029,-0.0151,-0.0319\right)}^T,\mathbf{\Psi }={\left({11.39}^{°},{27.60}^{°},{25.23}^{°}\right)}^T$.

So, $$\begin{array}{}\text{(52)}& {}_B{}^{I}g=\left[\begin{array}{cccc}0.8017& -0.3351& 0.4950& -0.3029\\ 0.3777& 0.9258& 0.0149& -0.0151\\ -0.4633& 0.1750& 0.8688& -0.0319\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

According to equation (11), there is $$\begin{array}{}\text{(53)}& {}_E{}^{I}g={}_B{}^{I}g{}_E{}^{B}g=\left[\begin{array}{cccc}0.7022& 0.3806& -0.6018& 3.5401\\ 0.2539& 0.6557& 0.7110& 0.2885\\ 0.6652& -0.6521& 0.3638& 1.5552\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

Situation 2.

Only the orientation of the base is obtained: $\mathbf{\Psi }={\left(11.39°,27.60°,25.23°\right)}^T$.

So, $$\begin{array}{}\text{(54)}& {{}^{I}R}_0=\left[\begin{array}{ccc}0.8017& -0.3351& 0.4950\\ 0.3777& 0.9258& 0.0149\\ -0.4633& 0.1750& 0.8688\end{array}\right].\end{array}$$

According to equation (14), there is $$\begin{array}{}\text{(55)}& {}_E{}^{I}g=\left[\begin{array}{cc}{{}^{I}R}_0{{}^{B}R}_E& {{}^{I}R}_0{{}^B\mathbf{p}}_E-\frac{{{}^{I}R}_0}{M}\sum _{i=1}^3{m}_i{{}^B\mathbf{p}}_i\\ \mathbf{0}& 1\end{array}\right]=\left[\begin{array}{cccc}0.7022& 0.3806& -0.6018& 3.5401\\ 0.2539& 0.6557& 0.7110& 0.2885\\ 0.6652& -0.6521& 0.3638& 1.5552\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

Situation 3.

Only the position of the base is obtained: ${{}^I\mathbf{p}}_0={\left(-0.3029,-0.0151,-0.0319\right)}^T$. $$\begin{array}{}\text{(56)}& {\mathbf{p}}_t=\frac{1}{M}\left(m_1{{}^B\mathbf{p}}_1,+m_2{{}^B\mathbf{p}}_2+m_3{{}^B\mathbf{p}}_3\right)={\left(0.2338,-0.0820,0.1779\right)}^T.\end{array}$$

According to equation (16), there is $$\begin{array}{}\text{(57)}& \left[\begin{array}{ccc}\text{c}\beta \text{c}\gamma & \text{s}\delta \text{s}\beta \text{c}\gamma -\text{c}\delta \text{s}\gamma & \text{c}\delta \text{s}\beta \text{c}\gamma +\text{s}\delta \text{s}\gamma \\ \text{c}\beta \text{s}\gamma & \text{s}\delta \text{s}\beta \text{s}\gamma +\text{c}\delta \text{c}\gamma & \text{c}\delta \text{s}\beta \text{s}\gamma -\text{s}\delta \text{c}\gamma \\ -\text{s}\beta & \text{s}\delta \text{c}\beta & \text{c}\delta \text{c}\beta \end{array}\right]\left[\begin{array}{c}0.2338\\ -0.0820\\ 0.1779\end{array}\right]+\left[\begin{array}{c}-0.3029\\ -0.0151\\ -0.0319\end{array}\right]=\mathbf{0}.\end{array}$$

So, $\mathbf{\Psi }={\left(\delta ,\beta ,\gamma \right)}^T={\left(11.385°,27.595°,25.236°\right)}^T$. Then, $$\begin{array}{}\text{(58)}& {}_B{}^{I}g=\left[\begin{array}{cccc}0.8017& -0.3352& 0.4949& -0.3029\\ 0.3778& 0.9257& 0.0150& -0.0151\\ -0.4632& 0.1749& 0.8688& -0.0319\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

According to equation (11), there is $$\begin{array}{}\text{(59)}& {}_E{}^{I}g={}_B{}^{I}g{}_E{}^{B}g=\left[\begin{array}{cccc}0.7021& 0.3805& -0.6019& 3.5399\\ 0.2541& 0.6558& 0.7109& 0.2891\\ 0.6652& -0.6520& 0.3637& 1.5556\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

4.2.2. Inverse Kinematical Solution of the Spatial Model

Likewise, final pose of EE and the iteration parameters are given at first. $$\begin{array}{}\text{(60)}& {}_E{}^{I}g=\left[\begin{array}{cccc}0.7022& 0.3806& -0.6018& 3.5401\\ 0.2539& 0.6557& 0.7110& 0.2885\\ 0.6652& -0.6521& 0.3638& 1.5552\\ & \mathbf{0}& & 1\end{array}\right].\end{array}$$

The max number of iterations is $N=100$. Convergence tolerance is $\epsilon =0.001$. Known parameters are $\mathbf{\Psi }={\left({11.39}^{°},{27.60}^{°},{25.23}^{°}\right)}^T$ (for Situation 2). ${{}^I\mathbf{p}}_0={\left(-0.3029,-0.0151,-0.0319\right)}^T$ (for Situation 3).

The initial joint angles are ${\mathbf{\theta }}_0={\left({45}^{°},{45}^{°},{45}^{°}\right)}^T$ and ${\mathbf{\Psi }}_0={\left(10,10,{10}^{\circ }\right)}^T$ (for Situation 3).

The results are listed in Table 3. Accuracy of the proposed iterative methods is demonstrated while applying to the tridimensional model.

Table 3

Simulation results of the Adams simulation and the iterations for the spatial model.

Method	Joint angles	Base’s attitude angles	Max error
Adams	$${\left(-{67.51}^{°},{88.35}^{°},{21.62}^{°}\right)}^T$$	$${\left({11.39}^{°},{27.60}^{°},{25.23}^{°}\right)}^T$$	—
Iteration for Situation 2	$${\left(-{67.504}^{°},{88.332}^{°},{21.645}^{°}\right)}^T$$	—	0.12%
Iteration for Situation 3	$${\left(-{67.263}^{°},{88.665}^{°},{21.874}^{°}\right)}^T$$	$${\left({11.590}^{°},{27.577}^{°},{25.065}^{°}\right)}^T$$	1.73%

5. Conclusions

This paper mainly studied the kinematical problem, including forward and inverse kinematics, of the free-floating space-robot system at the position level. Kinematical equations are systematically established based on the screw theory. Three situations with different known variables are discussed in detail. In particular, a numerical method for the inverse kinematical solution is proposed by using Newton’s iteration method. Examples demonstrate the accuracy and efficiency of the developed method, which appears to indeed be concise and valid compared to the traditional D-H method. The proposed iterative method can rapidly calculate joint angles with sufficient accuracy; it could be quite useful when differential kinematics is not concerned.

Conflicts of Interest

The authors declare that there is no conflict of interest regarding the publication of this paper.

Acknowledgments

This work is funded by the Shanghai Aerospace Science and Technology Innovation Foundation, Grant No. SAST2017-020.

References

[1] J. Denavit, R. Hartenberg, "A kinematic notation for low pair mechanisms based on matrices," ASME Journal of Applied Mechanics, vol. 22, pp. 215-221, 1955.

[2] A. Flores-Abad, O. Ma, K. Pham, S. Ulrich, "A review of space robotics technologies for on-orbit servicing," Progress in Aerospace Sciences, vol. 68,DOI: 10.1016/j.paerosci.2014.03.002, 2014.

[3] Z. Vafa, S. Dubowsky, "The kinematics and dynamics of space manipulators: the virtual manipulator approach," The International Journal of Robotics Research, vol. 9 no. 4,DOI: 10.1177/027836499000900401, 1990.

[4] W.-F. Xu, X.-Q. Wang, Q. Xue, B. Liang, "Study on trajectory planning of dual-arm space robot keeping the base stabilized," Acta Automatica Sinica, vol. 39 no. 1, pp. 69-80, DOI: 10.1016/S1874-1029(13)60008-7, 2013.

[5] Y. Umetani, K. Yoshida, "Continuous path control of space manipulators mounted on OMV," Acta Astronautica, vol. 15 no. 12, pp. 981-986, DOI: 10.1016/0094-5765(87)90022-1, 1987.

[6] K. Yoshida, Y. Umetani, "Control of space manipulators with generalized Jacobian matrix," Space Robotics: Dynamics and Control, pp. 165-204, DOI: 10.1007/978-1-4615-3588-1_7, 1993.

[7] W. Xu, Z. Hu, Y. Zhang, B. Liang, "On-orbit identifying the inertia parameters of space robotic systems using simple equivalent dynamics," Acta Astronautica, vol. 132, pp. 131-142, DOI: 10.1016/j.actaastro.2016.12.031, 2017.

[8] L. Yan, H. Yuan, W. Xu, Z. Hu, B. Liang, "Generalized relative Jacobian matrix of space robot for dual-arm coordinated capture," Journal of Guidance Control and Dynamics, vol. 41 no. 5, pp. 1202-1208, DOI: 10.2514/1.G003237, 2018.

[9] B. Liang, Y. Xu, M. Bergerman, "Mapping a space manipulator to a dynamically equivalent manipulator," Journal of Dynamic Systems, Measurement, and Control, vol. 120 no. 1,DOI: 10.1115/1.2801316, 1998.

[10] S. A. A. Moosavian, E. Papadopoulos, "On the kinematics of multiple manipulator space free-flyers and their computation," Journal of Robotic Systems, vol. 15 no. 4, pp. 207-216, DOI: 10.1002/(SICI)1097-4563(199804)15:4<207::AID-ROB3>3.0.CO;2-T, 1998.

[11] Y. Masutani, F. Miyazaki, S. Arimoto, "Sensory feedback control for space manipulators," Proceedings, 1989 International Conference on Robotics and Automation, pp. 1346-1351, DOI: 10.1109/robot.1989.100167, .

[12] C. Zhou, M. H. Jin, Y. C. Liu, Z. Zhang, Y. Liu, H. Liu, "Singularity robust path planning for real time base attitude adjustment of free-floating space robot," International Journal of Automation and Computing, vol. 14 no. 2, pp. 169-178, DOI: 10.1007/s11633-017-1055-1, 2017.

[13] X. Zhang, J. Liu, "Effective motion planning strategy for space robot capturing targets under consideration of the berth position," Acta Astronautica, vol. 148, pp. 403-416, DOI: 10.1016/j.actaastro.2018.04.029, 2018.

[14] B. Liang, W. F. Xu, Space Robotics: Modeling, Planning and Control, 2017.

[15] R. M. Murray, S. S. Sastry, Z. X. Li, A Mathematical Introduction to Robotic Manipulation, 1994.

[16] E. Sariyildiz, H. Temeltas, "A comparison study of three screw theory based kinematic solution methods for the industrial robot manipulators," 2011 IEEE International Conference on Mechatronics and Automation, pp. 52-57, DOI: 10.1109/icma.2011.5985630, .

[17] J. Gallardo, J. M. Rico, A. Frisoli, D. Checcacci, M. Bergamasco, "Dynamics of parallel manipulators by means of screw theory," Mechanism and Machine Theory, vol. 38 no. 11, pp. 1113-1131, DOI: 10.1016/S0094-114X(03)00054-5, 2003.

[18] K. Ibrahim, A. A. Ramadan, M. Fanni, Y. Kobayashi, A. A. A. Ismail, M. G. Fujie, "Kinematic analysis and control of limited 4-DOF parallel manipulator based on screw theory," 2012 IEEE International Conference on Systems, Man, and Cybernetics (SMC), pp. 47-52, DOI: 10.1109/ICSMC.2012.6377675, .

[19] C. R. Rocha, C. P. Tonetto, A. Dias, "A comparison between the Denavit–Hartenberg and the screw-based methods used in kinematic modeling of robot manipulators," Robotics and Computer-Integrated Manufacturing, vol. 27 no. 4, pp. 723-728, DOI: 10.1016/j.rcim.2010.12.009, 2011.

[20] F. Liu, H. Wu, "Research on modeling of robot manipulator dynamics based on screw theory," Journal of Jiangsu University Science and Technology (Natural Science Edition), vol. 22 no. 2, pp. 52-55, 2008.

[21] Y. P. Liu, H. T. Wu, "Deduction of Jacobian matrix for free-floating branching space robot based on recursion of screw," China Mechanical Engineering, vol. 21 no. 11, pp. 1275-1277, 2010.

[22] X. M. Zhang, Z. Z. Zhang, Y. B. Wang, T. Yang, T. Deng, "Inverse kinematical numerical method based on 6-DOF space manipulator," Missiles and Space Vehicles, vol. 3, pp. 81-84, 2016.