Pose Stabilization Control for Base of Combined System Using Feedforward Compensation PD Control During Target Satellite Transposition

TL;DR

This paper introduces a new control method to stabilize the base of a satellite transposition system using a dual-arm space robot.

Contribution

A novel pose stabilization method using feedforward compensation PD control for combined space systems is proposed.

Findings

The proposed method improves position and attitude accuracy compared to traditional control algorithms.

Simulation results show significant reduction in base pose errors under high-load and disturbed conditions.

Abstract

During the transposition of a target satellite, dynamic coupling between the target satellite, the manipulators, and the base frequently leads to disturbances in the base’s attitude. To deal with the issue, this paper proposed a pose stabilization method for the base of the post-capture combined system using the feedforward compensation PD control. Firstly, the mission sequence for repositioning a target satellite using a discrete-serpentine heterogeneous dual-arm space robot (DSHDASR) was analyzed. The dynamics model of the combined system, composed of the DSHDASR and a target satellite, was established based on the Newton–Euler recursive formulation. Then, the pose stabilization method integrating dynamic feedforward compensation and PD control was developed to stabilize the base of the combined system. Finally, the mission of target satellite transposition was simulated through the co-simulation model. Compared with the traditional control algorithms, the position accuracy and attitude accuracy for the proposed method showed an overall improvement. The results demonstrated that the proposed method significantly reduced base pose errors under high-load and disturbed conditions.

Figures

Full text

1. Introduction

Space robotic manipulators are essential for on-orbit servicing (OOS) [1,2], including debris removal, structure assembly, and satellite repair [3]. Based on the number of joints and link characteristics, space manipulators can be categorized into three types: discrete arms [4], serpentine arms [5], and continuum arms [6]. Unlike continuum arms, the discrete and serpentine manipulators present fewer challenges in terms of accurate modeling and precise control. Therefore, they are most widely adopted to carry out OOS. Furthermore, compared to a single discrete or serpentine arm, the DSHDASR combines the advantages of both types of arms, which ultimately overcomes specific limitations of individual manipulators [7]. Before conducting OOS, the space robotic arms installed on the DSHDASR should grasp the target satellite firstly, and then transpose it from the current location to the desired position for the subsequent operations [8]. Once the target satellite has been captured, the combined system is formed. However, during the target satellite transposition, strong dynamic coupling between the target, the manipulators, and the base can disturb the base’s pose (position and attitude). These disturbances not only reduce the tracking accuracy of manipulators but also bring potential risks to subsequent operations [9]. To validate space robotic control strategies under near-realistic microgravity conditions, various ground-based microgravity emulation platforms have been developed in recent years [10,11]. Therefore, suppressing these disturbances poses a great challenge for achieving precise tracking during target transposition, particularly under the condition of strong dynamic coupling.

To address this challenge, considerable research focused on pose stabilization control. Papadopoulos and Dubowsky [12] derived the dynamics model of space robotic manipulators using a Lagrangian framework, enabling coordinated control between the manipulator and the base. To achieve the desired attitude transposition, Yamada [13] proposed a joint path-planning method to achieve a desired change in attitude. By combining the Backstepping method and the radial basis function neural networks, Zhan et al. [14] designed an adaptive controller to ensure the trajectory tracking and base stabilization. Taking the end-effector as a Virtual Base, Zong et al. [15] established a new dynamics model of free-floating space manipulators to accomplish end-effector trajectory tracking and minimum disturbances to the base. Subsequently, Gong et al. [16] developed a dynamics model with parametric uncertainties to deal with the reactionless control of pose tracking. Based on the extended Jacobian model, Xu [17] designed a passivity adaptive controller to achieve the pose attitude adjustment and trajectory tracking. Considering the rigid body modes, Kuck et al. [18] developed a benchmark feedback controller to minimize trajectory tracking errors. While aforementioned methods can achieve base pose control, they do not account for the presence of a payload on the arm’s end-effector.

A payload significantly intensifies the dynamic coupling between the manipulators and the base. This brings a significant challenge to base stabilization control. To deal with the problem, Raina et al. [19] proposed a unified framework for impact dynamics modeling and post-capture reactionless control to guarantee attitude stabilization of the base. Wang et al. [20] developed a coordinated control method for post-capture operations, enabling a redundant space robot to maintain the attitude stability of the combined system. It should be noted that in both the methods above, the payload-to-base mass ratio is under 10%. For high proportion payloads, Chang et al. [21] presented an adaptive Backstepping-based nonsingular fast integral terminal sliding mode control method for the attitude stabilization of space robots after capturing a tumbling target. Woodward et al. [22] designed a Lyapunov-based control law to drive the base’s angular velocity attitude error asymptotically to zero. Based on dynamic coupling control with an iterative extended Kalman filter, Zhang et al. [23] proposed a rapid attitude stabilization strategy for highly robust base-attitude stabilization. Huang et al. [24] developed a robust attitude control method based on a disturbance observer and dynamic control allocation for the post-capture attitude control of the combined system, which ensured the prescribed attitude performance under uncertainties and external disturbances. Song et al. [25] applied a model-free control method for the attitude and orbital control of a post-capture combined spacecraft. Although these methods address base stabilization with large payloads on the end-effector, they do not consider the end-effector’s trajectory tracking performance.

Recently, the rapid development of intelligent algorithms has led to their application in base stabilization control. Liu et al. [26] investigated a trajectory planning and coordination control method to minimize the target’s attitude and disturbances to the base using the multi-objective particle swarm optimization (MOPSO) algorithm. Based on the deep Koopman operator-based modeling approach with model predictive control, Reference [27] achieved high-precision joint tracking and effective suppression of base attitude disturbances. According to reinforcement learning, Wang et al. [28] and Dong et al. [29] proposed a multi-target trajectory optimization method and a sequential optimization approach, respectively, to control the base’s dynamic response and stabilize its attitude. In addition, Ma et al. [30] applied deep reinforcement learning to restore a satellite’s attitude, which is a task analogous to base attitude control.

As discussed above, although numerous methods for base pose stabilization have been proposed, each method has certain limitations. Some approaches established the relationship between base pose disturbance and manipulator motion to achieve stabilization control of the base, yet did not account for the case of a payload on the end-effector. Other methods addressed base stabilization with a payload, but they did not account for end-effector pose tracking. Furthermore, while some intelligent algorithms could effectively handle base pose stabilization, their computational demands often made them impractical for real-time control on spacecraft hardware. To tackle these issues, a pose stabilization method based on feedforward compensation PD control was proposed for the base of the combined system, while the DSHDASR was transposing the target satellite from the initial pose to the desired pose. Compared with the PD control method and Backstepping method, this method introduces the reaction force and torque acting on the base in the control law, effectively decreasing the strong dynamic coupling effects on the base. In addition, compared with existing approaches based on adaptive schemes or intelligent algorithms, the proposed method has advantages including a clear computational structure, low computational complexity, and high computational efficiency. This feature makes the proposed strategy more suitable for spacecraft platforms within limited onboard computational resources and offers favorable potential for real-time implementation.

The remainder of this paper is organized as follows. Section 2 analyzes the overall process of the target satellite transposition. In Section 3, the dynamics model of the combined system is established by the recursive Newton–Euler formulation. Based on the PD control with feedforward compensation, Section 4 proposes the pose stabilization method for the base during the target satellite transposition. In Section 5, the proposed method is verified through the comparison simulations including the proposed method, the traditional PD method, and the Backstepping method. A comprehensive discussion of the challenges and simplifications is given in Section 6. Based on the analysis of previous sections, the last section summarizes the work of this study.

2. Overall Task Description

After capture, the target satellite often must be moved to a suitable location for subsequent on-orbit operations. This movement process is known as the target satellite transposition, as shown in Figure 1. At the initial state, the target satellite has been captured by the DSHDASR. The DSHDASR is composed of a base, a non-offset discrete arm (namely Arm-a) and a segmented cable-driven hyper-redundant serpentine arm (namely Arm-b). Together, the DSHDASR and the target satellite form the combined system. During transposition, Arm-a serves as the primary manipulator for moving the target satellite from the initial pose to the final pose (i.e., desired position and attitude) along a predefined trajectory, while Arm-b is used as the auxiliary manipulator, providing visual monitoring of Arm-a’s operation. Due to dynamic coupling between the base and the two manipulators, the base moves as Arm-a and Arm-b execute the coordinated transposition. As a result, the motion of the base alters the end-effector trajectories of both arms.

To ensure that the end-effectors of the two arms follow their predefined trajectories, it is necessary to maintain the base in a stable state. Therefore, this paper investigates control for stabilizing the base of the combined system during target satellite transposition. Some assumptions are given as follows.

• (1)The mass and geometric parameters of the target satellite and the DSHDASR are assumed to be known and are considered to be exact values in this study.
• (2)The base pose control can be achieved by thrusters and reaction wheels, whose installation positions and directions are known.
• (3)The target satellite and the end-effector of the non-offset discrete arm remain in a securely locked state during the transposing process.
• (4)Joint friction is neglected in the dynamic modeling and control design, and all joints are assumed to be ideal frictionless revolute joints.
• (5)The structural flexibility of the links and joints for the DSHDASR, as well as mechanical backlash and other non-ideal effects, are not taken into account.

3. Dynamics Modeling of the Combined System

As shown in Figure 2, the combined system is formed after the target satellite is captured by the DSHDASR.

As discussed above, the DSHDASR consists of a base, a non-offset discrete arm (Arm-a), and a segmented cable-driven serpentine arm (Arm-b). For Arm-a, each joint only has one rotational axis that is driven by a motor. The total number of joints is denoted by $[eqn]$ . In contrast, Arm-b is composed of $[eqn]$ segments with identical mechanical structures. Each segment has m associated universal joints, which are driven by several cables. An associated universal joint possesses a Pitch-axis and a Yaw-axis. For each segment, all associated universal joints are hybrid active–passive joints [31]. Specifically, the first universal joint in a segment is active and follows planned values, while the remaining joints are passive and move synchronously with it. {xI, yI, zI} and {x0, y0, z0} are the inertial frame and the frame fixed on the base center of mass, respectively. Unless otherwise stated, all the variables are represented in the inertial frame.

The Newton–Euler formulation is a well-known method that offers advantages in computational efficiency and numerical stability, making it suitable for real-time control. As a result, the Newton–Euler formulation is adopted to establish the dynamic model of the combined system.

3.1. Forward Recursion for Velocities and Accelerations Computation

As shown in Figure 2, $[eqn]$ and $[eqn]$ represent the ith joint and link of Arm-a. $[eqn]$ is the center of mass of $[eqn]$ and its position vector is denoted by $[eqn]$ . $[eqn]$ represents the vector from the joint $[eqn]$ to the center $[eqn]$ and $[eqn]$ is the vector from the center $[eqn]$ to the joint $[eqn]$ . For the ith joint of Arm-a, the position vector of $[eqn]$ is obtained as follows:

$$[eqn]$$

where $[eqn]$ is the position vector of the center of mass for the base. $[eqn]$ denotes the vector from the base’s center of mass to the first joint of Arm-a. According to the principle of relative motion, the angular velocity and the linear velocity of $[eqn]$ expressed in the ith joint frame are obtained as follows:

$$[eqn]$$

where $[eqn]$ stands for the transformation matrix from the (i − 1)th joint frame to the ith joint frame for Arm-a. Regarding the other parameters in Equation (2), the left superscript means that the vector is described in the frame indicated by the superscript, i.e., $[eqn]$ , $[eqn]$ , $[eqn]$ , and $[eqn]$ are all expressed in the (i − 1)th joint frame of Arm-a. The subsequent notations follow the same convention. $[eqn]$ and $[eqn]$ refer to the rotation axis and the angular velocity magnitude of the ith joint for Arm-a, respectively.

Differentiating Equation (2), the angular acceleration and linear acceleration of $[eqn]$ are derived as follows:

$$[eqn]$$

where $[eqn]$ denotes the angular acceleration magnitude of $[eqn]$ .

Based on Equation (3), the angular and linear accelerations of $[eqn]$ are obtained as follows.

$$[eqn]$$

For Arm-b, $[eqn]$ and $[eqn]$ are the jth joint and the jth link in the ith segment. $[eqn]$ refers the center of mass of $[eqn]$ and its position vector is named by $[eqn]$ . The position vector of the jth joint in the ith segment of Arm-b can be calculated as follows:

$$[eqn]$$

where $[eqn]$ is the vector from the base’s center of mass to the first joint of Arm-b. $[eqn]$ represents the vector from the joint $[eqn]$ to the center $[eqn]$ . $[eqn]$ is the vector from the center $[eqn]$ to the joint $[eqn]$ .

For the jth joint within the ith segment of Arm-b, its two rotational axes are denoted as $[eqn]$ and $[eqn]$ , and the corresponding angular velocities are $[eqn]$ and $[eqn]$ . Since the universal joints in the same segment are hybrid active and passive, the movement law of universal joints in the ith segment can be expressed as follows.

$$[eqn]$$

Based on the kinematics of relative motion, the angular velocity and the linear velocity of the jth joint in the ith segment are obtained as follows:

$$[eqn]$$

where $[eqn]$ is the transformation matrix from the second rotational axis of the (j − 1)th joint to the jth joint in the ith segment. With respect to the remaining parameters, the left superscript represents that the vector is expressed in the frame attached to the second rotational axis of the (j − 1)th joint in the ith segment. $[eqn]$ and $[eqn]$ are the first and second rotational axes of the jth joint. Their corresponding angular velocities are $[eqn]$ and $[eqn]$ , respectively.

Differentiating Equation (7), the angular and linear accelerations of $[eqn]$ are computed as follows.

$$[eqn]$$

According to Equation (8), the angular acceleration and linear acceleration of $[eqn]$ can be obtained as follows.

$$[eqn]$$

3.2. Backward Recursion for Joint Forces and Torques Calculation

The base’s motion including angular velocity, angular acceleration, linear velocity, and linear acceleration can be obtained from measurement by gyroscopes, accelerometers, or other sensors. In other words, $[eqn]$ , $[eqn]$ , $[eqn]$ , and $[eqn]$ are known. On the basis of Equations (4) and (9), the angular and linear accelerations of each link’s center of mass for Arm-a and Arm-b can be determined.

For both arms, joint forces and torques are calculated in reverse order, propagating from the end-effector back toward the base. For Arm-a, the dynamic relationship of the ith link is shown in Figure 3.

The gravitational force acting on each link of Arm-a can be considered as zero, due to the microgravity environment where the arm is. For the ith link of Arm-a, the force and torque balance equations are derived as follows.

$$[eqn]$$

where $[eqn]$ and $[eqn]$ are the inertia force and torque of the ith link. $[eqn]$ and $[eqn]$ refer to the force and torque exerted by the (i − 1)th link on the ith link. $[eqn]$ and $[eqn]$ are the force and torque exerted by the ith link on the (i + 1)th link. $[eqn]$ is the position vector from the ith joint to the (i + 1)th joint. They are all expressed in the frame fixed on the ith joint.

According to Newton’s laws and Euler’s equations of motion, the inertia force and torque of the ith link can be obtained as follows:

$$[eqn]$$

where $[eqn]$ and $[eqn]$ are the mass and inertia tensor of the ith link.

Substituting Equation (11) into Equation (10), the relationship between the motion and the applied force and torque for the ith link can be obtained as follows.

$$[eqn]$$

Here, $[eqn]$ and $[eqn]$ have been calculated based on Equation (4). Substituting Equation (4) into Equation (12), the following equations hold.

$$[eqn]$$

Notably, the issue of the virtual link $[eqn]$ has occurred during the backward recursion for calculating forces and torques. To address this, the external environment of Arm-a is treated as the link $[eqn]$ . If there is no contact between the end-effector and the external environment, the external force and torque at the end-effector are set to zero, i.e., $[eqn]$ and $[eqn]$ .

The recursive process for calculating the applied forces and torques on each joint of Arm-b is given as follows. The dynamic relationship of the jth link in the ith segment is shown in Figure 4.

The force and torque balance equations of the jth link in the ith segment are obtained as follows:

$$[eqn]$$

where $[eqn]$ and $[eqn]$ denote the force and torque exerted by the (j − 1)th link on the jth link. $[eqn]$ and $[eqn]$ are the force and torque exerted by the jth link on the (j + 1)th link. $[eqn]$ is the position vector from the jth joint to the (j + 1)th joint. They are all expressed in the frame attached to the second rotation axis of the jth joint.

Substituting Equation (9) into Equation (14), Equation (14) can be rewritten as follows.

$$[eqn]$$

Similarly, a virtual link appears in the above iteration, which does not exist for the actual structure of Arm-b. The treatment for this case is consistent with that of Arm-a.

4. Pose Stabilization of the Base by the PD Control with Feedforward Compensation

4.1. The Strategy for Pose Stabilization During Transposition

During target satellite transposition, the pose stabilization strategy for the base of the combined system is depicted in Figure 5. In this process, the non-offset discrete arm (Arm-a) is the mission manipulator to move the target satellite from its initial pose to its desired pose. At the same time, the segmented cable-driven hyper-redundant serpentine arm (Arm-b) serves as the auxiliary manipulator for monitoring the operation of Arm-a. Firstly, the desired pose of the target satellite was determined and the time (i.e., tf) was set based on transposition requirements. Accordingly, the target’s movement trajectory is planned using an interpolation method with a sampling time dt. After that, the pose of the target satellite is known at each moment.

Then, the system checked whether the current time exceeds the allotted mission time. If it is true, it means that the target satellite has been moved to the desired pose. Otherwise, the expected poses of end-effectors are calculated according to the geometric constraints between the target satellite and the two arms. Specifically, the relative poses between the target satellite’s center of mass, the capture point (Arm-a’s end-effector), and the monitoring point (Arm-b’s end-effector) remains fixed during transposition. Joint angular velocities and accelerations are determined by the inverse kinematics solution. The angular velocity and acceleration are obtained by the difference method. Furthermore, the PD control law integrated with dynamic feedforward compensation is designed to maintain the base stability. The feedforward compensation term is used to compensate for the base reaction force and torque generated by the motion of Arm-a and Arm-b, which can be calculated by the Newton–Euler recursion discussed in Section 3. The proportional-derivative control component is obtained according to the pose deviation of the base. Finally, the desired control force and torque are achieved by thrusters and reaction wheels installed on the base. This ensures the base remains stable while the DSHDASR performed the cooperative transposition-and-monitoring task.

4.2. Trajectory Planning of the Target Satellite

The schematic diagram for the target satellite transposition process is given in Figure 6. Before the transposition, the target satellite has been captured by Arm-a and the locked point is the capture point that is on the docking ring. During transposition, a monitoring device on Arm-b’s end-effector provided auxiliary monitoring of Arm-a’s operation. Therefore, the poses of Arm-a tip and Arm-b tip are dependent on the pose of the target satellite’s center of mass. Since the transposition is to move the target satellite from the initial pose to the desired pose, it can be taken as point-to-point motion planning. To ensure smooth start and stop transitions, the trajectory for the target satellite’s center of mass was generated using quintic polynomial interpolation.

The initial pose and the desired pose of the target are set as follows:

$$[eqn]$$

where $[eqn]$ , $[eqn]$ and $[eqn]$ are initial positions of the target satellite’s centroid along x-axis, y-axis, and z-axis, respectively. The desired positions are $[eqn]$ , $[eqn]$ , and $[eqn]$ . The initial attitude angles (Z-Y-X Euler angles) are denoted by $[eqn]$ , $[eqn]$ , and $[eqn]$ . Their desired attitude angles are $[eqn]$ , $[eqn]$ , and $[eqn]$ , respectively. The velocity and acceleration are set to zeros at the initial and terminal states. Therefore, the desired motion trajectory of the target satellite is obtained as follows.

$$[eqn]$$

According to the pose relationship among the target satellite’s centroid, Arm-a tip and Arm-b tip, the desired trajectories of the two arms are also calculated. Specifically, the desired poses of Arm-a tip and Arm-b tip are obtained, i.e., $[eqn]$ , $[eqn]$ , $[eqn]$ , and $[eqn]$ .

4.3. Feedforward Compensation PD Control Algorithm

As shown in Figure 7, the control method for the base’s pose stabilization is composed of two components. The former is dynamic feedforward compensation, while the latter is proportional-derivative control. The detailed calculation process for each component is given as follows.

4.3.1. Calculation of Feedforward Compensation

As discussed above, once the motion trajectory for the target satellite’s centroid was planned, the required end-effector poses for Arm-a and Arm-b became known. The joint angles, velocities, and accelerations of Arm-a and Arm-b are obtained by the inverse kinematics solution and numerical differentiation. These joint motion parameters are then fed into the Newton–Euler recursive algorithm to compute the corresponding joint forces and torques. Notably, the last link of Arm-a and the target satellite form a single rigid body after the target is captured. Therefore, the actual mass properties of the link used in the Newton–Euler recursion are changed, and they can be obtained as follows.

As shown in Figure 8, { $[eqn]$ , $[eqn]$ , $[eqn]$ } is the frame fixed on the last joint of Arm-a. { $[eqn]$ , $[eqn]$ , $[eqn]$ } is the frame attached to the target satellite. The equivalent centroid of the last link is denoted by $[eqn]$ . Its associated frame is { $[eqn]$ , $[eqn]$ , $[eqn]$ }. The position vector of $[eqn]$ expressed in { $[eqn]$ , $[eqn]$ , $[eqn]$ } is obtained as follows:

$$[eqn]$$

where $[eqn]$ and $[eqn]$ are the masses of the last link and the target, respectively. $[eqn]$ is the position vector from Arm-a tip to the centroid of the target satellite.

The equivalent inertia tensor expressed in the { $[eqn]$ , $[eqn]$ , $[eqn]$ } is obtained as follows:

$$[eqn]$$

where $[eqn]$ and $[eqn]$ are the inertia tensors of the last link and the target, respectively, expressed in their attached coordinate systems. $[eqn]$ is the unit matrix. $[eqn]$ is the rotation matrix from frame { $[eqn]$ , $[eqn]$ , $[eqn]$ } to frame { $[eqn]$ , $[eqn]$ , $[eqn]$ }. $[eqn]$ is the position vector from the centroid of the target satellite to $[eqn]$ .

Replace $[eqn]$ , $[eqn]$ , and $[eqn]$ with $[eqn]$ , $[eqn]$ , and $[eqn]$ . These equivalent parameters, along with the joint motion data, are then substituted into (13) and (15) to yield the joint forces and torques. Subsequently, reaction forces ( $[eqn]$ and $[eqn]$ ) and torques ( $[eqn]$ and $[eqn]$ ) on the base induced by the two arms are obtained. By expressing these forces and torques in the inertial frame, the dynamic feedforward compensation term is derived as follows:

$$[eqn]$$

where $[eqn]$ is the rotation matrix from the frame fixed on the first joint of Arm-a to the inertial frame. $[eqn]$ refers to the rotation matrix from the frame attached to the second rotational axis of the first joint of Arm-b to the inertial frame.

4.3.2. Control Law Design

As shown in Figure 7, the pose stability control of the base consists of the dynamic feedforward compensation term and the PD control term. The former is calculated using Equation (20). The PD control term is derived from the deviation between the base’s desired and measured poses. In practice, the desired pose is the base’s state when the target satellite is grasped. The base’s current pose is measured by onboard sensors (e.g., gyroscopes and accelerometers). Therefore, the PD control term is derived as follows:

$$[eqn]$$

where $[eqn]$ and $[eqn]$ are the position and attitude deviations for the base. Their linear and angular deviations are denoted by $[eqn]$ and $[eqn]$ . $[eqn]$ and $[eqn]$ are the proportional gain matrices corresponding to the position deviation and attitude deviation, respectively. $[eqn]$ and $[eqn]$ are the derivative gain matrices.

In this paper, the PD gains are designed following a systematic performance and robustness trade-off principle. Firstly, the gains should be proportional to the mass/inertia of the base, ensuring the physical consistency of the control force/torque. Furthermore, the gain matrix should be positive and definite to ensure system stability. Finally, to ensure system robustness, the minimum eigenvalue of the proportional gain matrix should be set in consideration of the magnitude of external disturbances and the required control accuracy. In other words, the proportional gain is expressed in a form explicitly related to the base mass, i.e., $[eqn]$ . Here, $[eqn]$ denotes the base mass, and $[eqn]$ represents a mass-independent normalized proportional gain that characterizes the desired closed-loop dynamic behavior.

Combining Equations (20) and (21), the control law of the base’s pose stabilization is obtained as follows.

$$[eqn]$$

4.3.3. Generation of Control Force and Torque

The control forces and torques calculated from Equation (22) are generated by thrusters and reaction wheels. As shown in Figure 9, the layout of the thrusters (see Figure 9a) and reaction wheels (see Figure 9b) installed on the base is illustrated. According to Figure 9a, twelve thrusters are mounted on the base, organized into three groups aligned with the x0, y0, and z0 axes of the base-fixed frame. Group 1 includes four thrusters that are numbered 1, 2, 3, and 4, respectively. These thrusters are arranged on the x0-z0 plane. Thrusters 1 and 2 produce thrust along the positive x0 axis direction, while Thrusters 3 and 4 generate thrust along the negative x0 axis direction. Group 2 consists of Thrusters 5, 6, 7, and 8, which are placed on the x0-y0 plane and Group 3 includes Thrusters 9, 10, 11, and 12, which are arranged on the y0-z0 plane. The definition of the thrust direction for these thrusters in Group 2 and Group 3 is the same as that of the thrusters in Group 1.

The position vector of the ith thruster in the frame fixed on the base centroid is denoted as $[eqn]$ . To generate a desired control force along a principal axis, two opposing thrusters from the corresponding group should be fired. For example, Thrusters 1 and 4 are selected to generate the desired control force along the x0 axis; Thrusters 5 and 8 are used for y0 axis; and Thrusters 9 and 12 are selected for the z0 axis. Thus, the desired control force from Equation (22) can be implemented as follows:

$$[eqn]$$ $$[eqn]$$ $$[eqn]$$

where $[eqn]$ is the thrust generated by the ith thruster. $[eqn]$ , $[eqn]$ , and $[eqn]$ are the components of the vector $[eqn]$ along the x0, y0, and z0 axes, respectively.

The layout of the reaction wheels is depicted in Figure 9b. Wheels 1, 2, and 3 are installed along the x0, y0, and z0 axes, respectively. Each wheel generates the control torque to adjust the attitude about its corresponding axis. Wheel 4 is oriented along an oblique axis to serve as a backup. In case of a failure in any one of Wheels 1–3, it is activated to maintain base attitude control.

Once the desired control torque in the feedforward compensation PD control algorithm is determined, it can be generated through rotation of the reaction wheels. Since the thrust direction of each thruster does not pass through the center of mass of the base, their operation introduces disturbance torques to the base. Therefore, the total torque output from the reaction wheels has two components, namely the nominal control torque from the algorithm and a compensatory torque to cancel thruster-induced disturbances. Taking Thrusters 1, 4, 5, 8, 9, and 12 as an example, the disturbance torque generated by their thrusts is obtained as follows.

$$[eqn]$$

Combining Equations (22) and (26) yields the total three-axis control torque required from the reaction wheels.

$$[eqn]$$

Assuming that Reaction Wheels 1, 2, and 3 are operating nominally, the relationship between the control torque and its angular acceleration is given as follows:

$$[eqn]$$

where $[eqn]$ , $[eqn]$ , and $[eqn]$ stand for moments of inertia of Reaction Wheels 1, 2, and 3, respectively. $[eqn]$ , $[eqn]$ , and $[eqn]$ refer to angular accelerations of Reaction Wheels 1, 2, and 3, respectively.

Substituting Equation (27) into Equation (28), the desired angular acceleration of each reaction wheel is calculated as follows.

$$[eqn]$$

Finally, the pose stabilization control of the combined system’s base is achieved by operating thrusters and reaction wheels according to Equations (23)–(25) and (29).

It is worth noting that actuator saturation can disrupt the controller’s linearity and may impair stability. To address this issue, the saturation limits are discussed as follows. From the perspective of control force and torque, their magnitude are determined by the feedforward compensation and the PD control term. Among these, the feedforward compensation contributes the majority, while the PD term accounts for a smaller portion. Therefore, to ensure that the feedforward compensation does not induce actuator saturation, a limiter has been applied to it. The specific settings are given as follows:

$$[eqn]$$ $$[eqn]$$

where $[eqn]$ and $[eqn]$ are the maximum output force and torque of the actuator.

4.4. Analysis of the Proposed Algorithm

4.4.1. Stability Analysis

The movement of the base is composed of the translation and rotation. According to Newton–Euler formulation, the dynamics model of the base can be obtained as follows:

$$[eqn]$$ $$[eqn]$$

where $[eqn]$ and $[eqn]$ are the mass and inertia of the base. Its linear velocity and angular velocity are $[eqn]$ and $[eqn]$ , respectively. $[eqn]$ and $[eqn]$ are the reaction force and reaction torque on the base induced by the two arms. $[eqn]$ and $[eqn]$ are the external disturbance force and torque, respectively.

Based on Equation (22), the force and torque of PD control with feedforward compensation are obtained. Substituting Equation (22) into Equations (32) and (33), it is calculated as follows:

$$[eqn]$$ $$[eqn]$$

where $[eqn]$ , $[eqn]$ . Because the base undergoes only small variations in its attitude, $[eqn]$ and $[eqn]$ can be calculated by $[eqn]$ and $[eqn]$ . $[eqn]$ and $[eqn]$ are the current position and attitude of the base. $[eqn]$ and $[eqn]$ are the desired position and attitude. Thus, $[eqn]$ and $[eqn]$ are the current linear velocity and angular velocity of the base. $[eqn]$ and $[eqn]$ are the desired linear velocity and angular velocity, and they are set as zeros in this paper. Therefore, $[eqn]$ , $[eqn]$ , $[eqn]$ , and $[eqn]$ can be obtained.

From a reformulation of Equations (34) and (35), the closed-loop dynamics for the base’s position and attitude can be derived as follows:

$$[eqn]$$ $$[eqn]$$

where $[eqn]$ and $[eqn]$ . $[eqn]$ and $[eqn]$ are defined as the lumped disturbance in the position channel and in attitude channel, respectively.

Assuming that $[eqn]$ and $[eqn]$ are zeros, Equations (36) and (37) can be rewritten as follows.

$$[eqn]$$ $$[eqn]$$

The Lyapunov function of the control system is established as follows.

$$[eqn]$$

Owing to the fact that $[eqn]$ , $[eqn]$ , $[eqn]$ , and $[eqn]$ are all positive definite matrix values, $[eqn]$ is greater than or equal to zero. Furthermore, only when $[eqn]$ , $[eqn]$ , $[eqn]$ , and $[eqn]$ are zeros, $[eqn]$ is equal to zero. Differentiating Equation (40), it can be obtained as follows.

$$[eqn]$$

Substituting Equations (38) and (39) into Equation (41), we can obtain

$$[eqn]$$

According to Equation (42), the control system is globally asymptotically stable.

The stability analysis presented above is derived under nominal conditions, where the dynamic feedforward compensation is assumed to exactly cancel the base coupling disturbances induced by the motion of the manipulators and the captured target satellite. Under this ideal assumption, the closed-loop system errors are guaranteed to asymptotically converge to zero.

In practical implementations, however, the perfect feedforward cancellation of coupling effects is difficult to achieve due to modeling uncertainties, parameter deviations, and unmodeled dynamics. In this case, the system is subject to bounded residual disturbances. Based on Lyapunov stability theory, it can be shown that the closed-loop system remains uniformly ultimately bounded, where the tracking errors converge to a bounded neighborhood whose size depends on the upper bound of the disturbances, rather than strictly converging to zero.

4.4.2. Robustness Analysis

Based on Equations (36) and (37), the lumped disturbances in position and attitude channels are $[eqn]$ and $[eqn]$ . They can both be decomposed into two components. The former is the feedforward compensation error, which could be generated due to parameter uncertainties or modeling errors. The latter is the disturbance force (or torque) exerted by the environment. When the base reaches the stable state, $[eqn]$ , $[eqn]$ , $[eqn]$ , and $[eqn]$ are all zeros. Therefore, its position and attitude errors can be obtained as follows.

$$[eqn]$$

Error bounds for the pose of the base are determined as follows:

$$[eqn]$$

where $[eqn]$ is the minimum eigenvalue of a matrix.

According to Equation (44), the position and attitude errors are related to external disturbances and proportional gain. Additionally, the larger minimum eigenvalue of the proportional gain results in a smaller error. As a result, we can enhance the robustness of the system by increasing the minimum eigenvalue of the proportional gain matrix.

4.4.3. Gain Margins Analysis

The system can be decoupled into six independent second-order systems.

$$[eqn]$$

For the position channel, the open-loop transfer function for each axis is given as follows.

$$[eqn]$$

The phase crossover frequency $[eqn]$ satisfies the following equation.

$$[eqn]$$

As derived from Equation (47), $[eqn]$ can be determined by using numerical methods such as the Newton–Raphson method or the bisection method. The corresponding gain at this frequency are calculated as follows.

$$[eqn]$$

Thus, the gain margin of the position channel is obtained as follows.

$$[eqn]$$

Similarly, the gain margin for the attitude channel can also be calculated, which is omitted here for brevity. It should be noted that the gain margin analysis presented in this section is conducted around the nominal operating point based on a linearized model. As such, it serves primarily as a theoretical stability reference and has inherent limitations when applied to the nonlinear and time-varying dynamics of the base–manipulator–target combined system, particularly under large payload conditions.

4.4.4. Sensitivity Analysis

During the transposition, the desired velocity and acceleration of the base are zeros, i.e., $[eqn]$ and $[eqn]$ . Substituting them into Equation (45), the closed-loop control system for the base is established as follows.

$$[eqn]$$

Based on Equation (50), the position or attitude channels can be treated as a second-order system that is fully defined by two characteristic parameters (damping ratio and natural frequency). Since the inherent parameters of a second-order system define its intrinsic dynamic characteristics, which remain unchanged by external disturbances, the disturbance force and torque can be neglected in Equation (50).

Taking the sensitivity analysis of parameter uncertainty as an example, assume that the true mass parameter of the base is $[eqn]$ . $[eqn]$ is the estimated mass of the base and $[eqn]$ is the uncertainty. From Equation (50), the closed-loop characteristic equation for the position channel can be derived as

$$[eqn]$$

From Equation (51), it can be seen that the uncertainty $[eqn]$ affects the distribution of the system poles. When $[eqn]$ increases (i.e., $[eqn]$ ), the poles move toward the imaginary axis, potentially making the response more oscillatory. In this case, the natural frequency decreases and the damping ratio reduces. When $[eqn]$ decreases (i.e., $[eqn]$ ), the poles move away from the imaginary axis, leading to an increase in both the natural frequency and damping ratio. However, if $[eqn]$ becomes too small, the poles might approach the imaginary axis again and exacerbate oscillations.

The damping ratio and natural frequency can be calculated as follows.

$$[eqn]$$

For the position channel, the sensitivity of the damping ratio and natural frequency with respect to the mass parameter uncertainty are then obtained as follows.

$$[eqn]$$

From Equation (53), it implies that a 10% increase in mass results in approximately a 5% decrease in both the damping ratio and natural frequency.

Similarly, the closed-loop characteristic function for the attitude channel can be established based on Equation (50). Its damping ratio and natural frequency can also be calculated. The corresponding sensitivity functions of the damping ratio and natural frequency with respect to the inertia parameter uncertainty are obtained as follows.

$$[eqn]$$

where $[eqn]$ and $[eqn]$ . A 10% increase in inertia reduces both the natural frequency and damping ratio by about 5%. This trend is fully analogous to the position channel analysis.

4.4.5. Discussion of the Challenges

Regarding the time delay issue, system delays from sensor sampling, communication, and computation could introduce phase errors, particularly degrading the performance of the model-dependent feedforward compensation. This would place a greater burden on the robustness of the PD feedback loop. We will study some strategies like predictive filtering or high-speed data protocols to address this issue.

On the problem of computational cost, the real-time execution of the recursive Newton–Euler algorithm for feedforward dynamics, while efficient (o(n)), must be rigorously validated on representative flight processor hardware.

About the actuator saturation, the physical limits of the thrusters and reaction wheels are the critical constraint. Each actuator requires integrating a saturation management block and a control allocation strategy to distribute torque among these actuators, preventing performance degradation or instability when commanded values exceed physical limits.

5. Simulation Study

To verify the proposed method, a co-simulation model is built using Matlab2016 and Adams2018 software. As shown in Figure 10, the DSHDASR for cooperatively transposing the target satellite is taken as an example in the simulation. At the initial moment of the simulation, the target has been grasped by the DSHDASR. The capture point is on the docking ring of the target, and it is referred to as the locked point. During the transposition, the discrete arm (Arm-a) is moving the target satellite from the initial pose to the desired pose, while the serpentine arm (Arm-b) is offering auxiliary monitoring for the operation of Arm-a.

5.1. The Parameters of Simulation Model

5.1.1. The Parameters of the DSHDASR

The DSHDASR is composed of the base, Arm-a, and Arm-b. Arm-a is a non-offset discrete manipulator with seven degrees of freedom (DoFs). Its D-H frames and parameters are given in Figure 11 and Table 1. The mass parameters are listed in Table 2.

Arm-b is a segmented serpentine arm with five identical segments. Each segment has six associated universal joints. All segments share an identical mechanical design. The D-H frames and parameters for a representative segment are shown in Figure 12 and provided in Table 3. The length and mass property of each link for Arm-b are given as follows.

$$[eqn]$$

The position and attitude (X-Y-Z Euler angle) of the base’s centroid with respect to the inertial frame are as follows.

$$[eqn]$$

The mass and inertia tensor matrix of the base are as follows.

$$[eqn]$$

The installation position and corresponding attitude (X-Y-Z Euler angle) for Arm-a and Arm-b, which are expressed in the frame attached to the base, are given as follows.

$$[eqn]$$

5.1.2. The Parameters of the Target Satellite

The mass and inertia tensor matrix of the target satellite are as follows.

$$[eqn]$$

When the target satellite is grasped, the position and attitude (X-Y-Z Euler angle) of the target satellite centroid with respect to the frame attached to the last joint of Arm-a are as follows.

$$[eqn]$$

5.2. The Desired Trajectory of the Target Satellite

To achieve target satellite transposition, the desired trajectory of the target satellite centroid is planned. Following Section 4.2, quintic polynomial interpolation ensured smooth acceleration and deceleration, generating the desired trajectory. In the simulation, the initial and final poses of the target satellite centroid are as follows.

$$[eqn]$$

Substituting Equation (61) into Equation (17), the variations of position and attitude for the centroid of the target satellite are obtained, as shown in Figure 13.

5.3. The Simulation for the Pose Stabilization of the Combined System Base

5.3.1. The Simulation Results for the Proposed Method

At the initial moment, the joint angles for Arm-a and Arm-b are set as follows.

$$[eqn]$$

From the target centroid trajectory (as indicated in Figure 13), the desired trajectories of the Arm-a tip and Arm-b tip can be obtained. Furthermore, the joint angles of both arms are calculated by the inverse kinematics solution. The variation curves of the joint angles are shown in Figure 14 and Figure 15, respectively.

In the simulation, the proportional and derivative coefficient matrices corresponding to the position and attitude of the base are set as follows.

$$[eqn]$$

Here, external disturbance forces and torques are introduced to simulate external environment influences, such as solar radiation pressure. The forces and torques are set as follows.

$$[eqn]$$

Figure 16 illustrates the base’s pose variations during transposition under the proposed control method. The maximum position deviations of the base centroid are 0.043 mm, −0.057 mm, and 0.06 mm along the xI, yI, and zI axes (see Figure 16a), respectively. The maximum deviations of attitude (X-Y-Z Euler angles, Figure 16b) are −0.006°, −0.024°, and 0.06°, respectively. The simulation results indicate that pose error is maintained within 0.1 mm and 0.1° by the proposed method, whose accuracy meets the engineering requirements.

Due to base pose deviations, the actual trajectory tracking of the target satellite could deviate from the desired trajectory. The trajectory tracking errors are shown in Figure 17. It can be observed that the maximum position errors of the target satellite centroid are −0.07 mm, 0.094 mm, and −0.094 mm (Figure 17a) corresponding to the xI, yI, and zI axes, while the maximum attitude errors during this process are 0.083°, –0.06°, and −0.072° (Figure 17b), respectively. Both position and attitude errors are limited within 0.1 mm and 0.1°, which further validate the correctness of the proposed method.

The control forces and torques, which are applied to the base, are shown in Figure 18. They can be generated by thrusters and reaction wheels. Figure 19 shows the 3D simulation models at the following milestones: t = 0 s, 12 s, 24 s, 36 s, 48 s, and 60 s.

5.3.2. The Simulation Results for the PD Control Method

To further validate the effectiveness of the proposed control strategy in stabilizing the base’s pose, a comparative simulation is conducted by using the PD control without feedforward compensation. When the traditional PD control method is applied, the pose variations for the base centroid of the combined system are given in Figure 20. As shown in Figure 20a, the simulation results indicate that the maximum position deviations of the base centroid along the xI, yI, and zI axes are 0.1289 mm, −0.1599 mm, and 0.1139 mm, respectively. It can be seen that all errors exceeded 0.1 mm. The maximum attitude deviations of the base centroid are −0.0186°, −0.0785°, and 0.1945° (see Figure 20b), respectively. Notably, the zI axis error exceeds 0.1°. In addition, the trajectory tracking deviations of the target satellite centroid are shown in Figure 21. It can be concluded that the maximum position tracking errors of the target satellite centroid are −0.2121 mm, 0.2665 mm, and −0.1857 mm, respectively. Meanwhile, they are 0.2442°, 0.1926°, and −0.2142° for the attitude tracking errors. All position errors and attitude errors are more than 0.1 mm and 0.1°, respectively. These results indirectly show the superior control accuracy of the proposed method over the traditional PD control method.

5.3.3. The Simulation Results for the Backstepping Method

In order to further verify the advantages of the proposed method, another comparison method (i.e., the Backstepping method) is implemented. The pose control errors of the base and the pose tracking error of the target satellite are shown in Figure 22 and Figure 23, respectively.

From Figure 22, it can be seen that the maximum position and attitude control errors of the base along xI, yI, and zI axes are 0.0636 mm, 0.1114 mm, 0.0724 mm, 0.0072°, 0.0381°, and 0.1141°, respectively. As given in Figure 23, the maximum position and attitude tracking errors of the target satellite along xI, yI, and zI axes are 0.1060 mm, 0.1840 mm, 0.1158 mm, 0.3176°, 0.1093°, and 0.2774°, respectively.

5.4. Simulation Results Analysis

The simulation results for the PD control method, the Backstepping method, and the proposed method are listed in Table 4.

According to the simulation results, the position control accuracy for the base of the combined system based on the proposed method is improved by 66.64%, 64.35%, and 47.32% along the xI, yI, and zI axes, respectively, compared with the PD control method. The attitude control accuracy of the base is improved by 66.13%, 69.55%, and 69.05%, respectively. Due to the dynamic coupling between the combined system’s base and the two manipulators, it will lead to the occurrence of a trajectory tracking error for the target satellite centroid. The proposed method achieves 67.04%, 64.62%, and 48.90% improvement in position tracking accuracy, compared with the PD control method. Additionally, they are 69.78%, 68.90%, and 69.84% for the attitude tracking accuracy.

Compared with the Backstepping method, the accuracy in both position and attitude control of the base are improved by 32.39%, 48.83%, 17.12%, 12.50%, 37.27%, and 47.24%, respectively. The position tracking accuracy of the target satellite are improved by 34.06%, 48.75%, and 18.05%, respectively. They are 76.76%, 45.20% and 76.71% for attitude tracking accuracy.

As a result, it can be concluded that the proposed method significantly improves the pose control precision of the base and the pose tracking accuracy of the target satellite. The effectiveness of the proposed method is thoroughly verified.

Regarding the comparison with the Backstepping method, it should be noted that the Backstepping controller is introduced in this study as a base stabilization strategy, with the same control objective as the conventional PD controller, namely suppressing the base disturbances induced by the manipulator and target satellite motion. The Backstepping method is not designed to directly optimize the attitude tracking performance of the target satellite. For a fair comparison, both controllers are tuned to achieve comparable settling times. In addition, relatively conservative virtual control gains are adopted for the Backstepping controller to ensure stability under strong dynamic coupling and large payload conditions.

In addition, the sensor noise is accounted for by incorporating model parameter uncertainty. The simulations are performed for parameter uncertainty levels of 2%, 5%, and 10%, respectively, to evaluate each scenario. The simulation results are listed in Table 5. It can be concluded that a larger parameter uncertainty leads to greater control errors. Furthermore, with a parameter uncertainty of 15%, the maximum position control error of the base is less than 0.09 mm, and the attitude error is below 0.11°. Meanwhile, for target satellite transposition, the position tracking error is less than 0.15 mm, and the attitude error remains under 0.12°.

In this analysis, the parametric uncertainty scenarios listed in Table 5 primarily focus on variations in the system mass and inertial properties. It is further clarified that sensor noise effects, which introduce measurement errors in the feedback loop and in the estimation of external disturbance forces and torques, are not explicitly modeled as stochastic processes, but are instead implicitly reflected through bounded uncertainty and equivalent disturbance terms within the adopted parametric uncertainty framework.

6. Discussion

A more comprehensive discussion of the challenges and simplifications associated with real-world implementation is given, including time delay, computational cost, and actuator saturation. Regarding the time delay issue, system delays from sensor sampling, communication, and computation could introduce phase errors, particularly degrading the performance of the model-dependent feedforward compensation. This would place a greater burden on the robustness of the PD feedback loop. We will study some strategies like predictive filtering or high-speed data protocols to address this issue. On the problem of computational cost, the real-time execution of the recursive Newton–Euler algorithm for feedforward dynamics, while efficient (o(n)), must be rigorously validated on representative flight processor hardware. About the actuator saturation, the physical limits of the thrusters and reaction wheels are a critical constraint. Each actuator requires integrating an actuator saturation management block and a control allocation strategy to dynamically distribute torque between actuators, preventing performance degradation or instability when commanded values exceed physical limits. It should be noted that the implementation-related aspects, including computation time per Newton–Euler recursion iteration, delay tolerance, and actuator saturation limits, are intended as qualitative considerations derived from simulation analysis. These discussions aim to provide preliminary insights into the practical feasibility of the proposed control strategy, rather than to report quantitatively validated hardware-level performance metrics. A detailed quantitative evaluation of these aspects requires hardware implementation or hardware-in-the-loop experiments, which are beyond the scope of the present study and will be addressed in future studies.

The base pose stabilization control method proposed in this study has so far been validated only in a simulation environment, and the analysis results are obtained under a set of idealized modeling assumptions. Specifically, the mass and geometric parameters of the system are assumed to be accurately known, while non-ideal factors such as joint friction, structural flexibility, and mechanical backlash are neglected in the simulation. Consequently, the simulation results mainly reflect the control performance of the proposed method in suppressing base dynamic coupling disturbances induced by large payload motions under ideal conditions. The advantage of simulation-based studies lies in their ability to evaluate the effectiveness and stability of the proposed strategy in a controllable and repeatable environment, which facilitates the investigation of system dynamics and the underlying control mechanisms. However, in practical physical implementations, factors such as sensor noise, actuator constraints, computational delays, and modeling inaccuracies may adversely affect control performance. To address these issues, future work will incorporate parameter uncertainties, joint friction, and structural flexibility, and will further assess the engineering applicability of the proposed method through higher-fidelity simulations or experimental platforms.

7. Conclusions

This paper proposed a pose stabilization algorithm based on feedforward compensation PD control to stabilize the base of the combined system during target satellite transposition. The overall transposition mission was analyzed to define the roles of each manipulator. Then, a dynamics model of the combined system was established to calculate the feedforward compensation term. A pose stabilization strategy combing feedforward compensation and PD control was developed to suppress disturbances to the base. Finally, to verify the proposed method, a co-simulation model was established byMatlab2016 and Adams2018 software. Comparative simulations were conducted for the proposed method and the traditional method. By incorporating the reaction forces and torques exerted on the base as feedforward compensation terms in the control law, the proposed method effectively mitigates the strong dynamic coupling effects during the target satellite transposition process. As a result, the errors of the base’s position and attitude are reduced significantly, and the trajectory tracking accuracy is improved for the transposition operation. Simulation results indicate that the proposed method exhibits good applicability for base stabilization control under large payload conditions, providing a feasible solution for maintaining base pose stability during the target satellite transposition process.

The current research only focuses on pose stabilization of the base for the combined system that is composed of a heterogeneous multi-arm space robot and a target satellite. The application of these methods for heterogeneous multi-arm systems has not been explored. Its dynamic recursive formulation does not account for joint friction or mechanical clearance. Furthermore, the robustness validation in this study is restricted to deterministic bounded-uncertainty scenarios rather than statistical characterization. Future work will incorporate stochastic analysis methods, such as Monte Carlo simulations, to statistically evaluate the control performance under random measurement noise and time-varying disturbances. In addition, the minimum energy optimization for target satellite transposition is also not considered in this paper. In the future, we will strive to do research on these topics.