Trajectory planning and base attitude restoration of dual-arm free-floating space robot by enhanced bidirectional approach

When free-floating space robots perform space tasks, the satellite base attitude is disturbed by the dynamic coupling. The disturbance of the base orientation may affect the communication between the space robot and the control center on earth. In this paper, the enhanced bidirectional approach is proposed to plan the manipulator trajectory and eliminate the final base attitude variation. A novel acceleration level state equation for the nonholonomic problem is proposed, and a new intermediate variable-based Lyapunov function is derived and solved for smooth joint trajectory and restorable base trajectories. In the method, the state equation is first proposed for dual-arm robots with and without end constraints, and the system stability is analyzed to obtain the system input. The input modification further increases the system stability and simplifies the calculation complexity. Simulations are carried out in the end, and the proposed method is validated in minimizing final base attitude change and trajectory smoothness. Moreover, the minute internal force during the coordinated operation and the considerable computing efficiency increases the feasibility of the method during space tasks.


Introduction
Extra vehicular activity plays an important role in space exploration and has been studied by many researchers. It mainly includes large space structure construction, satellite refurbishment and refueling, and space debris removal [1,2]. These tasks are dangerous and laborious for astronauts. As an alternative solution, the space robot attracts attention as a feasible, safe way for on-orbit tasks [3,4].
Among different varieties of space robots, the freefloating space robot (FFSR) stands out for the advantage of fuel saving [5,6]. During space tasks, the FFSR suffers from dynamic coupling that may result in satellite pose disturbance. The disturbance affects robot end positioning accuracy and may interfere with ground communication [7,8]. The end positioning problem can be well solved by the generalized Jacobian matrix-based control. Attentions are focused on minimizing satellite attitude disturbance because ground communication is remarkably influenced by satellite angular motion rather than translation motion [9].
To eliminate satellite base disturbance during trajectory planning, the reactionless null space method [10,11] and the null reaction force method [12,13] were exploited. The former plans trajectories with zero base angular velocity, whereas the latter minimizes the force exerted on the satellite base. Both methods are effective but can only be implemented for manipulators with sufficient redundant degrees of freedom (DOFs). The disturbance map method was proposed to plan low-disturbance paths [14], and the zero-reaction manipulator was designed to decouple the motion between the satellite and the manipulator [8]. However, these two methods cannot be expanded to high-DOF manipulators. The attitude control method was developed to control the base pose by adjusting the manipulators, but the method does not consider end tasks [15]. Free-falling cat was studied, and the near-optimal control approach was applied to find a minimum disturbance trajectory in joint space, but the method suffers from the velocity discontinuity and low computing efficiency [16]. Nakamura and Mukherjee [17] proposed the bidirectional approach (BA) to solve the nonholonomic problem. The method employs a virtual space robot as an assistant in trajectory planning and can reduce base attitude error while ensuring accurate joint positions. However, jumping joint velocities decrease the reliability and the feasibility.
In the aforementioned studies, researchers mainly focus on single-arm space robots that are fit for simple tasks. For complex tasks in space, multiarm space robots show better extensibility [18]. Multiarm space robots can accomplish tasks in more dexterous ways, such as coordinated operation, parallel task, and target clamping [1]. Studies on multiarm FFSR are also being carried out these years, including capturing tool designs [19], capturing strategies during coordinate operations [20], and satellite-manipulator decoupling strategies of freeflying and free-floating robots [1,21]. To eliminate base disturbance with the decoupling strategy, the rapidly exploring random trees-based method was applied [22]. However, this method needs additional DOFs and does not consider end tasks. Abdul Hafez et al. [1] optimized the satellite base movement of dual-arm robot with the quadratic optimization method. However, the velocity jumps conflict with the real manipulator control. The reactionless null space method was also applied to multiarm space robot to eliminate satellite base attitude disturbance [2,23]. However, additional DOFs or manipulators are needed to balance satellite base disturbance. The same problem was also met by Xie et al. [24], who tried to compensate the base disturbance with balance manipulator. These works mainly pay attention to dual-arm parallel tasks where the robot ends the movement in free space. Scenarios of coordinated operation with dual-arm FFSR are more complicated. Extra constraints are imposed to robot ends in coordinated operation, and the available DOFs are further reduced. As a result, more challenges are encountered for satellite base disturbance minimization.
The studies above show two ideas of solving the base-manipulator coupling problem: The first one tries to keep the base attitude stable in the whole dynamic process, which is mainly used in online trajectory planning, and the second one just focuses on the unity of the initial and final base attitude, which is mainly used in offline trajectory planning. Based on the second idea, this paper proposes a method named enhanced bidirectional approach (EBA) to solve the nonholonomic problem of dual-arm FFSR with limited redundant DOFs. The method aims to keep the satellite attitude the same at the initial and final moment when the dual-arm FFSR works in free space, and minimize the final attitude disturbance when the robot works with additional constraints. As an extension of BA, the method proposed in this paper aims to overcome velocity jumps in BA and ensure a velocitylevel continuous trajectory. In the EBA, the trajectory is planned in two situations: free end dual-arm robot and dual-arm robot with end constraints. The model of robot motion in free space is analyzed first, and a general state equation is then determined. A stable control law is introduced based on the analysis of the state equation, and the minimum disturbance trajectory is finally obtained. For tasks with end constraints, the constraint equations are integrated into the robot motion model, and the stable control law can be obtained directly in analogy with the situation of free end control. The EBA aims to decrease base disturbance effectively and ensures smooth joint trajectories. The smooth trajectory can avoid additional base attitude error and internal force in robot control.
The contents of this paper are organized as follows. The tasks and the FFSR analyzed are described in Section 2. The derivations of the EBA with and without end constraints are shown in Section 3. Simulations to validate the reliability and feasibility of the theory are carried out in Section 4. Finally, the paper is concluded in Section 5.
2 Problem description and robot model

Problem description
The dual-arm space robot is mainly composed of two space manipulators and a satellite base. When the space robot carries out on-orbit tasks, it usually runs into two scenarios, as shown in Fig. 1. In both figures, the opaque robot stays at the initial robot state, and the transparent pose is the final state. In Fig. 1(a), the two manipulators accomplish parallel tasks. Two manipulators are given separate, independent positioning tasks without contact, and the robot motion is planned in free space, similar to single-arm planning. In Fig. 1(b), the two manipulators accomplish tasks in a coordinated manner, such as capturing or recycling an object. In this situation, extra geometry and force constraints are imposed to manipulator ends; thus, further analysis is needed in motion planning.
The robot redundant DOFs are usually limited for satellite base disturbance elimination; thus, a real-time reactionless strategy is challenging. A reasonable solution is to plan the robot motion offline. In this manner, the concern is more about the initial and final states. Then, the planning goal becomes minimizing the final satellite attitude variation and ensuring precise positioning of the robot joints or ends.
In this paper, the task in free space is assumed to position the robot in the desired working pose, and the coordinated task is to capture and drag an object to the designated position. In both cases, the initial and final robot configurations (joint angles and satellite base attitude) are predefined, and our planning is to design smooth joint trajectories with eliminated final base attitude variation.

Robot model
The robot model and the corresponding D-H coordinate system are illustrated in Fig. 2. The robot contains a satellite base supporting two same-structured 7-DOF manipulators. The base coordinate system is located at the mass center of the satellite, and the installation coordinate systems of the manipulators are positioned at the opposite sides of the satellite. The D-H and installation parameters are listed in Table 1, and the inertia parameters are presented in Table 2. The mass centers of the links are measured with respect to the D-H coordinate systems, and the inertia is measured with respect to the corresponding mass centers.

Enhanced bidirectional approach
The EBA is proposed to find a smooth trajectory with restorable base attitude. In the method, the initial and final robot configurations are given first, and the trajectory is planned offline. To find a trajectory between the two configurations, a considerable idea is to seek trajectories from the two configurations simultaneously until the two trajectories meet each other. Suppose the trajectory starting from the initial state is T1 and the trajectory sought from the final one is T2. The planned trajectory can be defined as a combination T1 and the inverse of T2. In this manner, trajectory planning can be modelled as a convergence problem of the error between T1 and T2.
To propose a base attitude restoration method for dualarm FFSR, the dual-arm system is analyzed in two scenarios: dual-arm robot moving in free space and dualarm robot with end constraints. Moreover, the traditional BA method is introduced briefly in the section for a better  understanding of the EBA method.

Overview of the traditional bidirectional approach
Based on the said basic idea, trajectory planning can be equivalent to a "meeting motion" of the real robot and a virtual robot. The real robot is located at the initial configuration at the initial time as the opaque manipulator shown in Fig. 1, and the virtual robot stays at the final configuration at the initial time as the transparent one in Fig. 1. When the planning progress begins, the two robots start moving toward each other, and the process terminates when the state variables of the two robots coincide halfway. Then, T1 and T2 can be defined as the trajectories of the real and virtual robots, respectively.
To find the trajectory error between T1 and T2, the state equation is derived first. Based on the law of conservation of momentum, the satellite angular velocity of FFSR can be calculated as in Eq. (1): is the roll-pitch-yaw (RPY) angle of satellite base, α x , α y , and α z are respectively x, y, and z terms of the satellite base RPY angle, and are correspondingα α J sα J sω θ terms of and , respectively, and are the analytical Base-Jacobian and geometric Base-Jacobian, respectively, and is the joint angle.

α θθ
The robot configuration is defined as the state variable, that is, x = [ T , T ] T , and the input is defined by u = .
The state equation can be stated as follows: (2) where N is the number of robot joints, I N×N is the identity matrix with N dimension, and W is the input matrix of the robot system in the BA.
To calculate trajectory error, the state equations of the real and virtual robots are established based on Eq. (3). If the states of the real and virtual robots are noted by x 1 and x 2 , respectively, the following state equations can be obtained: {ẋ where W 1 and W 2 are the input matrices of the real and virtual robots in the BA, respectively, and u 1 and u 2 are the inputs of the real and virtual robots in the BA, respectively To find a convergent trajectory error, the state error Δx = x 1 -x 2 is selected to build the Lyapunov function, as shown in Eq. (4). The input ensuring a convergent trajectory error is derived as Eq. (5). With the input, the corresponding time-derivative of V is shown in Eq. (6).

ū W
where V is the Lyapunov function of the system, Δx = x 1 − x 2 is the state error between the real and virtual robots, Q is an arbitrary symmetric positive-definite matrix, is the augmented input composed by u 1 and u 2 , = [W 1 , −W 2 ], and (·) † means the pseudoinverse of corresponding matrix. Equations (4) and (6) show that the   V Lyapunov function V ≥ 0, the corresponding timederivative ≤ 0, and the equal sign holds if and only if Δx = 0. Thus, the system is asymptotically stable, and the steady-state value of V is zero.
The BA plans the trajectories of the real and virtual manipulators simultaneously. Supposing the two robot states coincide at meeting time t m , the real system input is combined as follows: W uW † ū W The derivation and simulation results of the BA show the following: 1) The desired joint angular velocities u(0) and u(2t m ) are nonzero, which results in velocity jumps at the initial and final time. 2) When state error Δx approaches zero, W 1 asymptotically equals W 2 . Thus, asymptotically tends to be singular. The simulation results show that though is the product of and Δx, in most cases, is divergent resulting from the singularity of . Then, the joint angular velocity jumps at time t m .
3) Owing to the velocity jumps above, when the BA is implemented in real robot control, the difference between the actual and desired trajectories results in extra base attitude variation.

Enhanced bidirectional approach for dual-arm robots without end constraints
The BA was proposed for the single-arm space robot. In practice, the BA can also be applied to dual-arm system without end constraints because both systems can plan the manipulator joint trajectories freely, that is, no additional constraint equations are needed if the BA is applied to the dual-arm system without end constraints. The EBA improves the BA by planning a smoother joint trajectory, and the situation of Fig. 1(a) is studied first. In the EBA, the system state equation is reconstructed, a new acceleration-level state equation is proposed first, and the corresponding system input is calculated and optimized for system stability and computing efficiency.

System model and input
The angular velocity discontinuities in the BA at time 0 (initial time), t m (meeting time), and 2t m (final time) suggest that only convergent Δx is not sufficient for successful robot control. Zero velocity is also needed at time 0, t m , and 2t m . Thus, an acceleration-level state equation is proposed in the EBA.
Joint angular acceleration is selected as input U, and then the new acceleration-level state equation for the FFSR can be defined as follows: x T ] T is the state variable of robot in the EBA, and z = is the joint angular velocity, x = [ T , T ] T represents a combination of the satellite RPY angle and robot joint angle, A and B are the state and input matrices of the system in the EBA, respectively, and W can be obtained by Eq. (2). Based on Eq. (8), the state equations of the real and virtual robots are listed as follows: where A 1 , A 2 and B 1 , B 2 are the state and input matrices of the real and virtual robots in the free-end system in the EBA, respectively, X 1 , X 2 and U 1 , U 2 are the state variables and the inputs of the real and virtual robots, respectively.

∈
To ensure the trajectory error between the real and virtual robots converging with zero joint velocities, trajectory error and joint velocity are combined to build the Lyapunov function. The combined variable s R 2N is introduced in Eq. (10), and the Lyapunov function is defined in Eq. (11).
∈ where P R 2N×(N+3) is an undetermined intermediate matrix that unifies the dimensions of Δx and z, and will be determined in the next section.
Taking a time derivative of Eq. (11) gives Equations (11) and (12) show that if the system input is calculated as in Eq. (13), is not larger than zero, as shown in Eq. (14).
U where k is an arbitrary positive number, and is the augmented input composed by U 1 and U 2 . The equal sign holds if and only if s = 0. Thus, the system is asymptotically stable, and the steady-state value of V is zero.
Substituting the joint angular velocity z into Eq.
is assumed, where m is an arbitrary positive number. Equation (16) is rewritten as ∆ẋ= − m∆x. (17) α θ α θẆ θW Define Δ as the base RPY angle error and Δ as the joint angle error, then Eq. (17) shows that Δx = [Δ T , Δ T ] T decays exponentially until around zero. If is a nonsingular matrix, z = will also converge to zero according to Eq. (15). However, similar to the BA, will gradually become singular as Δx approaches zero; thus, further input modification is required to keep the system stable.
where is the time derivative of and can be determined by Eq. (19).
is a high-order term of and is more susceptible to singularity [25].
W To eliminate the singularity influence on the system stability, is neglected in the paper, and the Moore-Penrose inverse is substituted by the damped least square inverse [26]. Then, Eq. (20) is obtained: z, (20) W † †W where is the damped least square inverse of and defined as . Then, the system input is derived as  The damped least square term may influence the system in two ways. A small damping factor λ has a minimalW † W † † influence on the system in the inception phase; thus, Δx still converges to zero. Though approaches infinity, still stays finite; hence, the first term in Eq. (20) converges to zero. The following Lyapunov function is considered: is positive definite, thus, V is asymptotic stable, and z gradually converges to zero.

α α θ
When damping factor λ is large, a large damping factor has a greater effect on the Base-Jacobian because W is mainly composed of the Base-Jacobian matrix and the identity matrix, where the former has a much lower order than the latter. The greater influence on the Base-Jacobian then leads to a fluctuated Δ . As a combination of Δ and Δ , Δx no longer approaches zero, and z fluctuates around 0 at time t m .
To sum up, a small damping factor causes overlarge joint angular acceleration and velocity, whereas a large damping factor results in large tracking errors. In practice, the damping factor should be selected according to the actual scenario.

Enhanced bidirectional approach for dual-arm robots with end constraints
The former section introduced the method for the singlearm or dual-arm space robot without end constraints. In this section, the trajectory planning problem in the dualarm coordinated operation is addressed. In this situation, the EBA is integrated with geometric constraints, and the corresponding input is modified. ω The dual-arm FFSR system is schematically illustrated in Fig. 3. O g is the inertial coordinate system; r a and r b are the end vectors of arms A and B, respectively; r ab is the vector pointing from arm B end to arm A end; v a , v b , a , Fig. 3 Structure of dual-arm free-floating space robot system. ω and b are the end velocities and the end angular velocities of arms A and B, respectively. All vectors are presented with respect to O g .
When the dual-arm robot is performing a coordinated operation, the velocities and angular velocities of the two arm ends satisfy the following constraints: { According to robot inverse kinematics, the constraints can be transformed into where H is the coefficient of the geometric constraints in coordinated operation, J Gva , J Gvb and , are the traditional end velocity and the angular velocity Jacobians of arms A and B, respectively, is the skew-symmetric matrix and is defined as θ Equation (26) shows that is located in the null space of H and where L is the null space of H, and ξ is an arbitrary vector. U = ξ is selected, and then Eq. (8) can be transformed into x where W L is the mapping matrix from the variable z to variable in the constraint-end robot system and is defined by W L = W•L, A L is the state matrices of constraint-end system in the EBA and B L is the input matrices of the constraint-end system in the EBA. Corresponding complete state equations are where A L1 and A L2 are the state matrices of the real and virtual robots in the constraint-end system in the EBA, and B L1 and B L2 are the input matrices of the real and virtual robots in the constraint-end system in the EBA. Then, in analogy with Eq. (20), the input can be derived as Eq. (31): where is the augmented mapping matrix in the constraint-end system, and can be calculated in analogy with Eq. (5).
As a null space of H, L is not rank full; thus, W L = W•L is a singular matrix. To limit the robot joint velocities, aW † † L large damping factor λ is needed in , and this factor selection results in the fluctuated Δα and z. To tackle this problem, a further modification is introduced to Eq. (31) to ensure the convergence of z: where Sig(t, t 0 ) is the Sigmoid function and can be obtained by Eq. (33), and t 0 is the time when velocities are desired to be zero. In practice, t 0 can be selected as the time when Δθ is located near zero.

Simulation and verification
The EBA aims to plan joint trajectories satisfying the joint angle demands and satellite base attitude demands. In this paper, two simulations are carried out: One aims to verify the control of the FFSR in the free space, and the other aims to validate the control of the dual-arm robot system with end constraints.
Each simulation is composed of two processes: planning and controlling. In planning, the initial and final conditions are first determined, and then the joint trajectories are planned. In controlling, the planned trajectories work as the joint controller inputs, and the joint space torque control method (closed-loop inverse dynamic control method for free space moving and master-slave for coordinated operation) is implemented for robot control. 4.1 Simulation of the free-floating space robot without end constraints The experimental and controlled groups are set in the simulation to validate the EBA methods. The BA, the 5degree-polynomial planning method, and the nearoptimal control approach are employed in the controlled group. The 5-degree-polynomial planning method is an offline method, and it plans the joint trajectory with 5degree-polynomial from the initial joint angle to the aim joint angle, as shown in Eq. (34). The near-optimal control approach is also an offline method, and it is proposed to optimize the final base attitude variation. In the method, the system state equation is also described by Eq. (2), and the system input is designed by a combination of finite number of Fourier orthogonal basis, as shown in Eq. (35). The coefficients of the orthogonal basis are solved by modified quasi-Newton algorithm during optimization [16]. In the controlling simulation, the traditional closed-loop proportional-differential (PD) inverse dynamic control method [27] is applied for robot control.
where t * is the initial time of trajectory planning, and the coefficients are calculated by the boundary condition of the trajectory.
where h is the number of the Fourier orthogonal basis, k cji and k sji are the coefficients of the ji terms of the Fourier orthogonal basis, and T is the total planning time of the near-optimal method.

Simulation parameter setting
The joint trajectory is planned according to the joint angles and the base attitudes at the initial and final times. Four methods are applied in the planning simulation, and the parameters are selected as follows. For the EBA, the parameters are set as k = 1.3, and m = 0.125. Given that k is larger than 10m, it is in accordance with k >> m and meets the demand of Eq. (20). Moreover, to eliminate the influence of the Sigmoid function and the damped least square term, t 0 = 1000 s and λ = 0 are set. For the BA, the effect of Q −1 on the system is similar to m in the EBA. Thus, Q = 8I is selected. For the near-optimal control approach, the number of Fourier orthogonal basis is 98 (=14×7), and the penalty factor is 10000.
The simulation time of planning is determined as follows. For BA and EBA, the time is determined by the value of Δx. The simulation time of 5-degree-polynomial planning method and the near-optimal control approach is affected by the joint velocity and joint acceleration. The unlimited-time simulation results show that the velocities of the BA are divergent, whereas the EBA works with a fine convergence and preferable stability. For an effective comparison, the simulation time of the planning simulation is set as 300 s. Within this period, meeting velocities of the EBA converge to 0.001°/s while the condition of the BA is acceptable. The simulation time of the two other controlled groups is set as 20 s, which limits the joint velocity and joint acceleration within 10°/s and 10°/s 2 .
In the controlling simulation, the traditional closed-loop PD inverse dynamic control method is applied. In the method, the PD parameters are selected according to the joint acceleration. For the feasibility of the control system, excessively large acceleration should be eliminated to enable motor driving. In the BA method and near-optimal control approach, small PD parameters are selected to limit joint acceleration because the velocities are discontinuous. In the EBA and the 5-degreepolynomial planning method, continuous velocities permit much larger PD parameters. In the simulation, the acceleration is limited by 10°/s 2 ; thus, PD parameters are selected as proportional parameter K p = 0.15, differential parameter K d = 0.6 for the velocity jumping group, and K p = 10, K d = 40 for the continuous velocity group. In the controlling simulation, the simulation time is extended to show the step response process at the velocity jumping moment.

Simulation results and analysis
Though the simulation time of the 5-degree-polynomial planning method and the near-optimal control approach is shorter than that of BA and EBA, the actual calculation time of BA and EBA is approximately 370 s, and the actual run time of the 5-degree-polynomial planning method and the near-optimal control approach are approximately 0.42 and 6300 s, respectively. The run time of the 5-degree-polynomial planning method is much less than that of the others because the calculation is achieved only by a polynomial computation. This process is fast, but the method does not optimize the base attitude. For the other methods, the BA and the EBA have higher calculation efficiency than the near-optimal control approach.
The joint angular velocities of simulations are illustrated in Fig. 4. Figures 4(a)-4(d) are the joint angular velocities of the EBA, the BA, the 5-degreepolynomial planning method, and the near-optimal control approach, respectively. To compare the BA and the EBA clearly, their joint angular velocities are dispersed along the axis labelled as "Robot joint order," as shown in Figs. 4(a) and 4(b). In both figures, joint orders 1-7 represents the joints 1-7 of arm A, and orders 8-14 mean the joints 1-7 of arm B. The controlling simulation time for the BA, the EBA, and the nearoptimal control approach is extended. For the BA and the near-optimal control approach, extensions are used to show the complete the step response process, and the extension of the EBA is used as a comparison of BA.
In Fig. 4, the desired joint angular velocity provided by the planning simulation is plotted by dashed lines, and the real joint angular velocity obtained from the controlling simulation is depicted by solid lines. Figure 4(a) shows that the real joint angular velocity of the EBA fits the desired values well. The velocities at the initial and final times of all joints are zero, and at the meeting time, the joint velocities are also near zero. Moreover, the "meeting points" are marked by filled circle in the figure, the curves before the points are trajectory T1, and the curves after the point are the trajectory T2. In Fig. 4(b), the desired joint angular velocity of BA jumps at t = 0 s, t = 150 s, and t = 300 s. These sudden changes occur at the initial, meeting, and final times, and induce step responses to real curves. For the 5-degree-polynomial planning method shown in Fig. 4(c), the real curves coincide with the desired curves well. As a well-known second-order continuous method, the controlling simulation time is not expanded, and the results prove that the robot will stop moving at the final time without step response. Figure 4(d) shows that the desired joint angular velocity is continuous during the planning except at the initial and final times. The discontinuity of the desired trajectory causes a step response and tracking error in the real trajectory and decreases the reliability of the method. In the figure, the real (solid lines) and desired curves (dashed lines) are plotted to show the complete responses in the real robot control (especially the step response for the discontinuous desired trajectory). The step response introduces tracking errors to the system; thus, the final base attitude deviates from its theoretical value.
The joint angles of simulations are illustrated in Fig. 5 to show the robot working state, and the real and desired trajectories are plotted for a complete response presentation. Figures 5(a)-5(d) are the joint angles of the EBA, the BA, the 5-degree-polynomial planning method, and the near-optimal control approach, respectively. As in Fig. 4(a) Though the final state error of the near-optimal control approach is slightly larger than those of the three other methods, it is still located in the neighbor of the required angle and satisfies the initial and final conditions. For BA, the curve peaks of Fig. 5(b) suggest that the method works with less controlling precision than the three other methods. The figures show that all the methods meet the general requirement of ensuring precise robot arm positioning. The advantage and disadvantage of the methods can be concluded from the trajectory continuity, final base attitude change, and computing efficiency.
The trajectories of satellite base attitudes are plotted in Fig. 6. Figures 6(a), 6(b), and 6(d) are the base RPY angles of the EBA, the BA, the 5-degree-polynomial planning method, and the near-optimal control approach, respectively.   the largest one because no optimization is applied in this method. This phenomenon also suggests the effectiveness of the other optimization methods. For the EBA, the real and desired final base attitudes are the smallest. It proves that the method can decrease the final base attitude variation effectively, and the smooth trajectory has a minimal perturbation on the system.

Simulation of free-floating space robot with end constraints
The experimental and controlled groups are set in the simulation to validate the EBA method. The last section proved that the trajectories planned by the BA perform worst in robot control. If the method is applied to robot tasks with end constraints, the situation is even worse. Thus, in the section, BA is not selected in the controlled group. To fit the coordinated operation, the leaderfollower method is applied instead of the 5-degreepolynomial planning method [28]. In the method, one arm (leader) is planned by the 5-degree polynomial as in Eq. (34), and the other one (follower) is planned by velocitylevel inverse kinematics with Jacobian matrix, that is, the controlled group consists of the leader-follower method and the near-optimal control approach. In the planning simulation, the simulation time of the EBA is 200 s, and the parameters are selected as k = 1.1, m = 0.1, t 0 = 95 s, and λ = 10 −5 . When t 0 is 95 s, the meeting angles converge to 0.01°, which is in accordance with the requirement in Eq. (32). For the near-optimal control approach and the leader-follower method, the simulation time is 20 s, the number of Fourier orthogonal basis is 98 (=14×7), and the penalty factor is 10000. As mentioned above, the simulation time and parameters are selected based on the joint velocity and the acceleration limits.
In the controlling simulation, additional position and force constraints are imposed. The position constraints are satisfied by the planning method, and the force constraints are satisfied by the master-slave method [29]. Given that the EBA plans the robot motion in the joint space, joint space impedance control is applied to the slave arm, and the impedance parameters are selected as M = 1000, B = 400000, and K = 10. Moreover, the PD parameters in the master-slave method affect the end positioning precision. To obtain a precise relative end pose, the PD parameters for the EBA and the leader-follower method are K p = 1000 and K d = 800. For the near-optimal control approach, as a result of velocity discontinuity, the parameters are selected as K p = 0.15 and K d = 0.6.

Simulation results and analysis
The actual run time of the leader-follower method, the EBA, and the near-optimal control approach are approximately 22, 270, and 6200 s, respectively. The planning efficiency of the EBA is lower than that of the leader-follower method but higher than that of the nearoptimal control approach.
Figures 7(a)-7(c) are the joint angular velocities of the EBA, the leader-follower method, and the near-optimal control approach, respectively. The real and desired curves are plotted to indicate the complete system responses. In Fig. 7(a), the real joint angular velocity of the EBA coincides with the desired one, and both start from zero and end at zero. Moreover, the transition velocity of the EBA is also zero, and the transition period is marked by the filled circles in the figure. For the leader-follower method, the initial and final velocities are also zeros. The joint velocities of the two arms show different performances because the leader arm and the follower arm are planned in different ways. In the nearoptimal control approach, the desired initial and final velocities are not zero, and step responses exist in the curves. As a result, the real velocity does not fit the desired one well, and the controlling simulation time is extended to 25 s to show the step response. In the coordinated operation, the step responses also introduce trajectory tracking error, and it not only results in additional variation to the base attitude control but also has an influence on the system internal force.
The joint angles of the simulations are shown in Fig. 8. . For the EBA and the near-optimal control approach, the initial and final angles are all located in the neighbor of the required values (though the error of the near-optimal control approach is slightly larger than that of the EBA). For the leader-follower method, because the joint motion of arm A is planned by the 5-degree polynomial and joint motion of arm B are planned by the Jacobian-based method, the final joint angle of arm A satisfies the initial and final conditions, whereas the joint angle of arm B does not.  However, the end pose of the two arms are the same as those of the other methods. The angles in Fig. 8 show that the robot ends are located at the neighbor of the desired pose relative to the base at the final time.
The satellite base RPY angles are plotted in Figs Further simulation results show that the base attitude of the EBA can be reduced more when a smaller damping factor λ is selected. However, the motion range would further expand and eventually exceed the defined joint angle scope with a smaller λ. The base attitude of the leader-follower method becomes largest because the method does not have the ability of optimization. Though the desired final base attitude of the near-optimal control approach is zero, the trajectory tracking error greatly increases the real value. As a result, the real final attitude of EBA is still the smallest one.
For the coordinated operation, a large internal force may damage the robot system; thus, the force control results are focused here as an important item under supervision to show the influence of tracking error on the system thoroughly. In this paper, the master-slave method is applied to decrease the robot internal force. In the method, the master arm is controlled by the rigid body  Fig. 10 Impedance force of slave arm in dual-arm system without end constraints: (a) impedance force in the enhanced bidirectional approach, (b) impedance force in the leader-follower method, and (c) impedance force in the near-optimal control approach. inverse dynamics, and the slave arm is controlled by the joint space impedance control method. The internal impedance torque of the slave arm is shown in Fig. 10. For the EBA, the internal torques are all below 0.6 N•m, and the values can be even limited within 0.1 N•m if the start and end moments are excluded. In Fig. 10(b), the internal joint impedance torques are slightly larger, and the upper limits of the torques are approximately 2 N•m. Simulations show that increased planning time decreases internal torque. For example, if planning time is 35 s, the upper limit of the internal torque will be 0.25 N•m. Compared with the two other methods, the internal joint impedance torques of the near-optimal control approach are much larger, and the upper limit is more than 10 N•m. The internal force of the near-optimal control approach is the largest one, and this phenomenon results from step responses in the tracking trajectories during the controlling.
In summary, the EBA is compared with the four methods in the whole simulations: the BA, the 5-degreepolynomial planning method, the leader-follower method, and the near-optimal control approach. The EBA shows the best performance in final base attitude optimization and joint trajectory tracking, and has the minimum internal force in the coordinated operation. The 5-degree-polynomial planning method and the leaderfollower method work worst because they plan the joint trajectory without optimization. The trajectory discontinuity of the BA and the near-optimal control approach introduce additional attitude variation and internal force, thus proving the necessity of the improvement of the BA on velocity-level continuous trajectory. Moreover, compared with the other methods, the EBA has a considerable computing efficiency and is suitable for base attitude restoration in space tasks.

Conclusions
To enable the FFSR to start and end with a desired configuration while restoring the base attitude, the EBA is investigated in this paper. By reformulating the joint trajectory planning into an acceleration-level "meeting problem", a robust system model is proposed. The EBA is an off-line planning method that aims to achieve robot free-space tasks and dual-arm robot coordinated operations. Simulations are carried out with and without end constraints. The results verify that compared with the other methods, the EBA decreases the final base variation effectively and ensures smooth joint trajectories when generating robot joint trajectories. The smooth trajectories avoid additional base attitude error and internal force in the robot control; thus, it ensures positioning accuracy and control stability, and benefits the robot in system control, end positioning, and communication. Moreover, the high computing efficiency makes the method practical in real space tasks.