Landing control method of a lightweight four-legged landing and walking robot

Abstract

The prober with an immovable lander and a movable rover is commonly used to explore the Moon’s surface. The rover can complete the detection on relatively flat terrain of the lunar surface well, but its detection efficiency on deep craters and mountains is relatively low due to the difficulties of reaching such places. A lightweight four-legged landing and walking robot called “FLLWR” is designed in this study. It can take off and land repeatedly between any two sites wherever on deep craters, mountains or other challenging landforms that are difficult to reach by direct ground movement. The robot integrates the functions of a lander and a rover, including folding, deploying, repetitive landing, and walking. A landing control method via compliance control is proposed to solve the critical problem of impact energy dissipation to realize buffer landing. Repetitive landing experiments on a five-degree-of-freedom lunar gravity testing platform are performed. Under the landing conditions with a vertical velocity of 2.1 m/s and a loading weight of 140 kg, the torque safety margin is 10.3% and 16.7%, and the height safety margin is 36.4% and 50.1% for the cases with or without an additional horizontal disturbance velocity of 0.4 m/s, respectively. The study provides a novel insight into the next-generation lunar exploration equipment.

Key words： landing and walking robot; lunar exploration; buffer landing; compliance control

1 Introduction

The probers, including lander and rover, play a pivotal role in scientific exploration on extraterrestrial bodies, such as Moon and Mars [1]. Luna-2, the first hard lander to reach lunar surface, can only perform a transient probe before colliding with the moon. Probers with landing capability are produced to realize surface exploration on the Moon; they usually adopt three or four foldable legs with invert-tripod mechanism, such as Surveyor-1 [2] and Luna-16 [3], or with cantilever mechanism, such as Apollo project [4] and Chang’e 3, 4, and 5 [5–7]. For probers with invert-tripod or cantilever mechanism, their leg often consists of a main branch chain and two multifunctional auxiliary branch chains. All branch chains use internal aluminum honeycomb material to absorb landing impact energy. However, the exploration capacity of such probers is still limited by their immovable nature. The drill or dig devices [8] mounted on landers can only sample the soil surrounding the fixed landing site, which may possibly be polluted by engine propellant. Therefore, researchers have proposed two main strategies to improve the exploration capacity of probers.

The first strategy is to use the cooperative work mode of lander and rover, such as the cooperation of Apollo lander and lunar roving vehicle (LRV) [9], and Chang’e lander and Yutu-1 and Yutu-2 [5,6]. However, all the existing rovers, such as LRV and Yutu-1 and Yutu-2, are wheeled mobile robots. Although they possess excellent mobility on relatively even terrain, their traversing performance on irregular terrains, such as gully, slope, and rock, remains poor even under active suspension system. At present, the legged robots, which are driven by hydraulic actuator or electric motor, attract increasing attentions owing to the excellent mobility performance on many challenging terrains [10]. Compared with the hydraulic actuator that is fit for legged robots working on the Earth, such as BigDog [11] and HyQ2Max [12], an electric actuator is considered to be more realistic to use on the Moon. For example, the hexapod SpaceClimber [13], whose legs driven by 24 motors were arranged symmetrically along the sagittal plane of the body, was designed to execute extraterrestrial exploration. This probe succeeded to provide mobility in slopes with inclinations of up to 40° and capacity of negotiating obstacles of up to 40 cm height. An eight-legged scorpion driven by 32 motors [14] was further designed by this group to enhance the terrain adaptability. Electric quadruped robots, such as MIT Cheetah 3 [15], ANYmal [16], and SpaceBok [17,18], are more widely studied, owing to their less driving modules. The SpaceBok was designed for roving on the Moon; it used parallel elastic elements to store and release energy and mounted a reaction wheel to adjust the pitch angle for achieving continuous jumping gait. These legged robots show good agility and traversability in complex terrain. However, the collision between body and ground is prone to occur during landing due to the short buffer stroke and low deceleration ratio in leg design, which will damage the recoil engine nozzle installed at the bottom of the robot.

Movable landers, which can perform a buffer landing and complete walking on the lunar surface, are developed [19]. A legged mobile lander using four identical 2-UPS&UP parallel mechanism was designed by Lin et al. [20], which integrated the functions of landing, walking, deploying, and orientation adjusting. To improve the landing impact performance of a movable lander, Zhou et al. [21] proposed a newly designed gradual energy-absorbing structure subjected to impact loads. Lin et al. [22,23] addressed a novel synthesis method for legged mobile landers using truss-mechanism transformation based on the structure of Chang’e lander. Han et al. [24,25] reported the dimensional design and synthesis methods for cantilever-type legged lander using performance atlas. All these mobile landers have not completed the prototype design, and their buffer landing performance needs to be further tested. The walking function of these mobile landers is achieved by adding an additional motor to each leg branch chain, with impact absorbed by internal aluminum honeycomb material during landing, which leaves the problem of irreversible deformation of the buffer material unsolved. Thus, these landers are still only capable of one-off buffer landing.

The above two traditional strategies to improve the probe performance still have some limitations. Scientists have proposed a novel method to achieve long-range mobility by flying on a planet with a thin atmosphere. The conceptual design of Dragonfly lander [26] was proposed to achieve landing on the Titan, in which four rotors were used and many instruments, including sampling drill, cameras, and multimission radioisotope thermoelectric, were carried. This conceptual design was proved to be feasible as the Ingenuity helicopter [27], with a weight of 1.8 kg, a height of 0.5 m, and a peak motor power of 350 W, successfully achieved multiple takeoffs and landings on Mars. This flying approach can be borrowed to the Moon exploration. As shown in Fig.1, the lander utilizes the recoil engines to hover in the air and then uses the vision system to find a relatively flat terrain for buffer landing. After completing the exploration task near the current landing site as a rover, the recoil engines will rework and drive the lander to fly over complex terrains, such as valleys or mountains, and the above procedures will be repeated to achieve buffer landing on another site. On the basis of this assumption, a repetitive hexapod robot integrated with walking and landing functions is designed by our team. However, this hexapod robot has more driving modules and larger real-time communication data with controller compared with the four-legged lander. These features require more powerful recoil engines and more fuel, resulting in a lethal problem of higher error rate, so the detailed landing ability in the cases of leg failure was systematically discussed in our previous study [28].

Fig.1 Conceptual sketch of repetitive landing.

Full size|PPT slide

In this study, a four-legged landing and walking robot called “FLLWR” is designed for lunar exploration. It integrates the functions of repetitive lander and rover on a lightweight body and has less communication data simultaneously. The impact energy is converted into heat energy on motors through hybrid compliance to avoid the damage of internal structure, in which a passive spring is mounted between the thigh and shank to provide passive compliance, and three integrated drive units (IDUs) are mounted on each leg to provide active compliance. During the landing period, each IDU imitates the dynamics of an active torsional spring and an active torsional damper by a compliance control method [29,30]. Our robot is promising to achieve repetitive flying and landing by recoil engines because this approach has protected the intact functions of legs, especially in challenging terrains, such as deep craters or mountaintops on the Moon [31]. The landing control method is a critical technical challenge for the FLLWR because it determines the success of landing. A novel method with landing method framework, state switching strategy, trajectory planning, compliance planning, and height and attitude stabilizer, is proposed to solve this problem. The results of numerical simulations and experiments on a five-degree-of-freedom lunar gravity testing platform (5-DoF-LGTP) without/with horizontal disturbance velocity simultaneously show that this novel method can be used to achieve buffer landing for the FLLWR.

The rest of this paper is organized as follows. Section 2 introduces the robot design. Section 3 establishes the kinematics and dynamics of quadruped robots. Section 4 proposes the landing control method. Section 5 gives the landing numerical simulation. Section 6 discusses the experimental results on the 5-DoF-LGTP. Section 7 provides the conclusions.

2 Four-legged landing and walking robot

As shown in Fig.2, an FLLWR is designed to implement three functions: actively controlled landing, roving exploration, and repetitive landing. Its size is about 1.35 m length, 1.35 m width, and 0.70 m height, and its weight is about 40 kg. The robot mainly consists of one body and four identical legs. All legs are evenly and symmetrically distributed under the body. The body mainly carries the control system, sensing system, power systems, and other payloads. The control system consists of motor drivers, I/O modules, and motion controller. The sensing system includes inertial measurement unit (IMU), laser radar, encoders, torque sensors, and others. The power system includes battery packages. Load weights are used to emulate other payloads for lunar exploration. The body is covered by a square pyramid shell made of aluminum alloy, and the internal electronics are protected by the shell.

Fig.2 Structure of a four-legged landing and walking robot.

Full size|PPT slide

2.1 Multifunction design

The robot has multiple functions, including folding, unfolding, repetitive landing, and walking. The demonstrations for these functions are given in Fig.3. The robot keeps folding morphology during the launching period on the Earth with all legs folded under the body, and this compact structure reduces the occupied volume in rocket. When the robot reaches the predetermined landing site, the propulsion engine begins to decelerate and stabilize the attitude of the robot until the landing conditions are met. When the robot is about to land on the Moon, the unfolding process is initiated to extend all legs for lengthening compression stroke and enhancing support stability. The propulsion engine is turned off, the robot begins free falling, and the landing control method starts work at the same time. Each leg is controlled in accordance with the foot-ground touch events and the foot-ground contact forces.

Fig.3 Multifunction design during different phase: (a) buffer landing like a cantilever lander and (b) roving like a four-legged robot.

Full size|PPT slide

The robot becomes a four-legged rover adopting a stable wave gait [32] with at most one leg swinging and at least three legs supporting after successful landing. The walking process for one gait period is shown in Fig.3(b). In locomotion, the green arrow line indicates the center of mass (COM) projection movement from the old green point to the new red point. The red arrow arc indicates that one leg swings, and the triangle formed by blue dotted line indicates the supporting polygon.

2.2 Leg design

The prototype and mechanism of the leg are presented in Fig.4. The leg mechanism is a hybrid parallel-serial mechanism, where a parallel five-bar mechanism (PFBM) is serially connected with a revolute pair (RP). The PFBM is a parallelogram planar parallel mechanism, where the connecting link is parallel with the thigh, and the rocker arm is parallel with the shank. Each leg has three active degrees-of-freedom (DoFs) and is actuated by three IDUs. The PFBM is used to realize the swing and stretching/shrinking movement of the leg in the leg sagittal plane. The RP between the body and the PFBM is used to realize the abduction/adduction movement of the leg in the body longitudinal plane. A passive linear spring is added between the thigh and the shank to balance the partial weight of the lander and reduce the required driving torque and the impact on robot from foot-ground event.

Fig.4 Three active degrees-of-freedom hybrid leg design: (a) leg prototype and (b) leg mechanism. IDU: integrated drive unit, RP: revolute pair, PFBM: parallel five-bar mechanism.

Full size|PPT slide

2.3 Actuator design

An external rotor motor rather than an internal one is selected to provide more powerful drive capability in the same volume. A reducer is required in the actuator to provide the high torque during landing. The motor and the reducer are then integrated and designed as an IDU, consisting of servo motor, reducer, torque sensor, encoder, brake, bearing, coupler, and shell. Fig.5 shows the prototype of IDU and the main components. In accordance with the torque and velocity requirements, the motor parameters are selected as follows: peak torque of 2 N·m and maximum speed of 5000 r/min, the reducer parameters are selected as transmission ratio of 100 and maximum output torque of 200 N·m. For the feedback control, the rotary encoder with 4096 counts per turn is selected. The torque sensor has a precision of 0.1% and a measurement range of ±250 N·m. The IDU is compact in structure and weighs 1.1 kg.

Fig.5 Design of the integrated drive unit: (a) three-dimensional and (b) inner structure.

Full size|PPT slide

2.4 Control system

The hardware structure of the control system is illustrated in Fig.6. The control algorithm runs on a real-time operating system in the controller. The roll and pitch angles of the body are sensed by IMU at 0.5 kHz. On the basis of the EtherCAT bus system at 4 kHz, the actual torque {\mathit{\tau }}_{\mathrm{a}\mathrm{j}}, actual position {\mathit{\theta }}_{\mathrm{a}\mathrm{j}}, and actual velocity {\dot{\mathit{\theta }}}_{\mathrm{a}\mathrm{j}} of each joint are detected by the torque sensor and encoder. Remote communication is used to achieve the human–robot interaction. The 12 command joint torques {\mathit{\tau }}_{\mathrm{c}\mathrm{j}} are transmitted to the IDUs by the EtherCAT bus as well.

Fig.6 Hardware structure of motion controller. RTOS: real time operating system.

Full size|PPT slide

3 Kinematics and dynamics of quadruped robots

3.1 Coordinate frames

Several coordinate frames are defined to establish the robot model, as shown in Fig.7. For the body coordinate frame {\mathrm{\Sigma }}_{\mathrm{b}}, the origin O_{\mathrm{b}} is located at the body center, x_{\mathrm{b}} is the forward direction, y_{\mathrm{b}} is the left direction, and {\mathit{\text{z}}}_{\mathrm{b}} is the upper direction. For the world coordinate frame {\mathrm{\Sigma }}_{\mathrm{w}}, its directions coincide with {\mathrm{\Sigma }}_{\mathrm{b}} at the initial moment, and the origin O_{\mathrm{w}} is located at the projection point of O_{\mathrm{b}} on the ground. For the ith leg coordinate frame {\mathrm{\Sigma }}_{\mathrm{l}i}, the origin O_{\mathrm{l}i} is located at the center point of the hip joint, x_{\mathrm{l}i} points to the body center, and y_{\mathrm{l}i} is parallel to {\mathit{\text{z}}}_{\mathrm{b}}. For the IMU coordinate frame {\mathrm{\Sigma }}_{\mathrm{i}\mathrm{m}\mathrm{u}}, the axis directions are the same as those of {\mathrm{\Sigma }}_{\mathrm{b}}.

Fig.7 Coordinate system in motion planning: (a) robot coordinate and (b) leg coordinate.

Full size|PPT slide

3.2 Leg modeling

The leg model includes leg forward kinematic model, leg inverse kinematic model, velocity model, and force model. The leg forward kinematic model is used to establish the relationship between the joint angles and the tiptoe position in the leg coordinate frame. The joint angles form a generalized coordinate vector \mathit{q}=[\begin{array}{ccc}{\theta }_{\mathrm{s}}& {\theta }_{\mathrm{t}}& {\theta }_{\mathrm{r}}\end{array}]{}^{\mathrm{T}}, where {\theta }_{\mathrm{s}} is the side angle, {\theta }_{\mathrm{t}} is the thigh angle, and {\theta }_{\mathrm{r}} is the rocker arm angle. The tiptoe position {\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}}=[\begin{array}{ccc}{x}_{\mathrm{t}\mathrm{i}\mathrm{p}}& y_{\mathrm{t}\mathrm{i}\mathrm{p}}& {\mathit{\text{z}}}_{\mathrm{t}\mathrm{i}\mathrm{p}}\end{array}{]}^{\mathrm{T}} can be obtained as follows:

(1){\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}}=\left[\begin{array}{c}{L}_{\mathrm{t}}\mathrm{c}\mathrm{o}\mathrm{s}{\theta }_{\mathrm{t}}-L_{\mathrm{s}}\mathrm{c}\mathrm{o}\mathrm{s}{\theta }_{\mathrm{r}}\\ \mathrm{c}\mathrm{o}\mathrm{s}{\theta }_{\mathrm{s}}\left(L_{\mathrm{t}}\mathrm{s}\mathrm{i}\mathrm{n}{\theta }_{\mathrm{t}}-L_{\mathrm{s}}\mathrm{s}\mathrm{i}\mathrm{n}{\theta }_{\mathrm{r}}\right)\\ \mathrm{s}\mathrm{i}\mathrm{n}{\theta }_{\mathrm{s}}\left(L_{\mathrm{t}}\mathrm{s}\mathrm{i}\mathrm{n}{\theta }_{\mathrm{t}}-L_{\mathrm{s}}\mathrm{s}\mathrm{i}\mathrm{n}{\theta }_{\mathrm{r}}\right)\end{array}\right],

where L_{\mathrm{t}} and L_{\mathrm{s}} are the lengths of thigh and shank, respectively.

The leg inverse kinematic model is employed to calculate the joint angles from the tiptoe position in the leg coordinate frame. For the given tiptoe position, we have

(2){\theta }_{\mathrm{s}}=\arctan \frac{{\mathit{\text{z}}}_{\mathrm{t}\mathrm{i}\mathrm{p}}}{y_{\mathrm{t}\mathrm{i}\mathrm{p}}},

(3){\theta }_{\mathrm{t}}=\arcsin \frac{-L_{\mathrm{s}}\sin B}{\sqrt{x_{\mathrm{t}\mathrm{i}\mathrm{p}}^2+y_{\mathrm{t}\mathrm{i}\mathrm{p}}^2/\phantom{y_{\mathrm{t}\mathrm{i}\mathrm{p}}^2{\cos }^2{\theta }_{\mathrm{s}}}{\cos }^2{\theta }_{\mathrm{s}}}}+\arctan \frac{y_{\mathrm{t}\mathrm{i}\mathrm{p}}}{x_{\mathrm{t}\mathrm{i}\mathrm{p}}\cos {\theta }_{\mathrm{s}}},

(4){\theta }_{\mathrm{r}}={\theta }_{\mathrm{t}}-B,

where

(5)B=\arccos \frac{-x_{\mathrm{t}\mathrm{i}\mathrm{p}}^2-y_{\mathrm{t}\mathrm{i}\mathrm{p}}^2/\phantom{-x_{\mathrm{t}\mathrm{i}\mathrm{p}}^2-y_{\mathrm{t}\mathrm{i}\mathrm{p}}^2{\cos }^2{\theta }_{\mathrm{s}}+L_{\mathrm{t}}^2+L_{\mathrm{s}}^2}{\cos }^2{\theta }_{\mathrm{s}}+L_{\mathrm{t}}^2+L_{\mathrm{s}}^2}{2L_{\mathrm{t}}{L}_{\mathrm{s}}}.

In accordance with Eqs. (2)‒(4), the generalized vector of the input angles \mathit{q}=[\begin{array}{ccc}{\theta }_{\mathrm{s}}& {\theta }_{\mathrm{t}}& {\theta }_{\mathrm{r}}\end{array}]{}^{\mathrm{T}} is obtained.

Let \dot{\mathit{q}}=[\begin{array}{ccc}{\dot{\theta }}_{\mathrm{s}}& {\dot{\theta }}_{\mathrm{t}}& {\dot{\theta }}_{\mathrm{r}}\end{array}{]}^{\mathrm{T}} be the generalized vector of the joint velocity. The tiptoe velocity {\dot{\mathit{P}}}_{\mathrm{t}\mathrm{i}\mathrm{p}}=[\begin{array}{ccc}{\dot{x}}_{\mathrm{t}\mathrm{i}\mathrm{p}}& {\dot{y}}_{\mathrm{t}\mathrm{i}\mathrm{p}}& {\dot{\mathit{\text{z}}}}_{\mathrm{t}\mathrm{i}\mathrm{p}}\end{array}]{}^{\mathrm{T}} by differentiating Eq. (1) can be expressed as

(6){\dot{\mathit{P}}}_{\mathrm{t}\mathrm{i}\mathrm{p}}={\mathit{J}}_{\mathrm{v}}(\mathit{q})\dot{\mathit{q}},

where {\mathit{J}}_{\mathrm{v}}(\mathit{q}) is the velocity Jacobian matrix as follows:

(7)\begin{array}{rl}& {\mathit{J}}_{\mathrm{v}}(\mathit{q})=\\ & \left[\begin{array}{ccc}0& -L_{\mathrm{t}}\sin {\theta }_{\mathrm{t}}& L_{\mathrm{s}}\sin {\theta }_{\mathrm{r}}\\ -\sin {\theta }_{\mathrm{s}}\left(L_{\mathrm{t}}\sin {\theta }_{\mathrm{t}}-L_{\mathrm{s}}\sin {\theta }_{\mathrm{r}}\right)& L_{\mathrm{t}}\cos {\theta }_{\mathrm{s}}\cos {\theta }_{\mathrm{t}}& -L_{\mathrm{s}}\cos {\theta }_{\mathrm{s}}\cos {\theta }_{\mathrm{r}}\\ \cos {\theta }_{\mathrm{s}}\left(L_{\mathrm{t}}\sin {\theta }_{\mathrm{t}}-L_{\mathrm{s}}\sin {\theta }_{\mathrm{r}}\right)& L_{\mathrm{t}}\sin {\theta }_{\mathrm{s}}\cos {\theta }_{\mathrm{t}}& -L_{\mathrm{s}}\sin {\theta }_{\mathrm{s}}\cos {\theta }_{\mathrm{r}}\end{array}\right].\end{array}

The length of passive spring can be written as follows:

(8)L_{\mathrm{p}\mathrm{s}}=\sqrt{L_{\mathrm{r}}^2+L_{\mathrm{f}}^2-2L_{\mathrm{l}}{L}_{\mathrm{r}}\cos (\text{π}-{\theta }_{\mathrm{t}}+{\theta }_{\mathrm{r}})},

where L_{\mathrm{l}} and L_{\mathrm{r}} are the distances from the intersection point of thigh and shank to the left and right fixed points of spring, respectively, and L_{\mathrm{p}\mathrm{s}} is the current length of the passive spring.

Its elastic potential energy of passive spring (E_{\mathrm{p}\mathrm{s}}) is

(9)E_{\mathrm{p}\mathrm{s}}=\frac{K_{\mathrm{p}\mathrm{s}}{\left(L_{\mathrm{p}\mathrm{s}}-L_{\mathrm{p}\mathrm{s}0}\right)}^2}{2},

where K_{\mathrm{p}\mathrm{s}} and L_{\mathrm{p}\mathrm{s}0} are the stiffness coefficient and original length of the passive spring, respectively.

The virtual work of passive spring (\text{δ}{W}_{\mathrm{p}\mathrm{s}}) is

(10)\text{δ}{W}_{\mathrm{p}\mathrm{s}}=-\frac{\mathrm{\partial }{E}_{\mathrm{p}\mathrm{s}}}{\mathrm{\partial }\mathit{q}}\text{δ}\mathit{q},

where \text{δ}\mathit{q} is the virtual displacement of the generalized coordinate vector \mathit{q}.

The relationship between the tiptoe force {\mathit{F}}_{\mathrm{t}\mathrm{i}\mathrm{p}}= [\begin{array}{ccc}{F}_x& F_y& F_{\mathit{\text{z}}}\end{array}{]}^{\mathrm{T}} and the joint torque vector \mathit{\tau }=[\begin{array}{ccc}{\tau }_{\mathrm{s}}& {\tau }_{\mathrm{t}}& {\tau }_{\mathrm{r}}\end{array}{]}^{\mathrm{T}} can be obtained by using the d’Alembert virtual work principle as follows:

(11){\mathit{F}}_{\mathrm{t}\mathrm{i}\mathrm{p}}^{\mathrm{T}}\text{δ}{\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}}=\left({\mathit{\tau }}^{\text{T}}-{\mathit{H}}^{\text{T}}\right)\text{δ}\mathit{q},

(12){\mathit{H}}^{\mathrm{T}}=[\begin{array}{ccc}0& h& -h\end{array}],

(13)h=\frac{-K_{\mathrm{p}\mathrm{s}}\left(L_{\mathrm{p}\mathrm{s}}-L_{\mathrm{p}\mathrm{s}0}\right)L_{\mathrm{l}}{L}_{\mathrm{r}}\sin ({\theta }_{\mathrm{t}}-{\theta }_{\mathrm{r}})}{L_{\mathrm{p}\mathrm{s}}},

where \text{δ}{\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}} is the virtual displacements of tiptoe and active joint, and H is the correction vector.

Considering the derivative relationship in Eq. (6), Eq. (11) can be rewritten as follows:

(14)\mathit{\tau }={\mathit{J}}_{\mathrm{v}}(\mathit{q}{)}^{\mathrm{T}}{\mathit{F}}_{\mathrm{t}\mathrm{i}\mathrm{p}}+\mathit{H}.

An additional correction vector H is added to the joint torque vector τ due to the use of passive spring. When the leg is compressed for buffering, the passive spring reduces the thigh torque amplitude and increases the shank torque amplitude. The peak torque is reduced, and the buffering capacity is improved because the amplitude of thigh torque is greater than that of shank torque in the buffering process.

3.3 Body modeling

The body model includes kinematic model and dynamic model. The kinematic model is used to express the relationship of position vector between the leg coordinate frame {\mathrm{\Sigma }}_{\mathrm{l}i} and the world coordinate frame {\mathrm{\Sigma }}_{\mathrm{w}}. As shown in Fig.7, {}^{\mathrm{l}i}{\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}} is the tiptoe position in the leg coordinate frame {\mathrm{\Sigma }}_{\mathrm{l}i}. {}^{\mathrm{b}}{\mathit{P}}_{\mathrm{l}i} is used to denote the origin position O_{\mathrm{l}i} of {\mathrm{\Sigma }}_{\mathrm{l}i} in the body coordinate frame {\mathrm{\Sigma }}_{\mathrm{b}}. {}^{\mathrm{b}}{\mathit{R}}_{\mathrm{l}i} is the rotation transformation matrix from {\mathrm{\Sigma }}_{\mathrm{l}i} to {\mathrm{\Sigma }}_{\mathrm{b}}. {}^{\mathrm{w}}{\mathit{P}}_{\mathrm{b}} is the body position in the world coordinate frame {\mathrm{\Sigma }}_{\mathrm{w}}. {}^{\mathrm{w}}{\mathit{R}}_{\mathrm{b}} is the rotation transformation matrix from {\mathrm{\Sigma }}_{\mathrm{b}} to {\mathrm{\Sigma }}_{\mathrm{w}}. The tiptoe position {}^{\mathrm{w}}{\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}} in the world coordinate frame is given as follows:

(15){}^{\mathrm{w}}{\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}}={}^{\mathrm{w}}{\mathit{P}}_{\mathrm{b}}+{}^{\mathrm{w}}{\mathit{R}}_{\mathrm{b}}\left({}^{\mathrm{b}}{\mathit{P}}_{\mathrm{l}i}+{}^{\mathrm{b}}{\mathit{R}}_{\mathrm{l}i}{}^{\mathrm{l}i}{\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}}\right).

A virtual dynamic model is established for the gravity compensation and the attitude adjustment. The mass of the leg can be ignored because it is smaller relative to the mass of the body to simplify the dynamics model. The contact between the tiptoe and the ground is inelastic, and no slipping or rebounding occurs. As shown in Fig.7, {\mathit{r}}_{\mathrm{c}\mathrm{o}\mathrm{m}} is the vector from the COM to the origin O_{\mathrm{b}}. {\mathit{r}}_i is the vector from the origin O_{\mathrm{b}} to the ith tiptoe. {\mathit{F}}_i is the ith tiptoe force from the ground. The body dynamics can be formulated as

(16)m_{\mathrm{b}}{\mathit{a}}_{\mathrm{b}}=m_{\mathrm{b}}\mathit{g}+\sum _{i=1}^4{\mathit{F}}_i,

(17)\frac{\mathrm{d}}{\mathrm{d}t}\left({\mathit{I}}_{\mathrm{b}}{\text{ω}}_{\mathrm{b}}\right)=\sum _{i=1}^4\left(\left({\mathit{r}}_{\mathrm{c}\mathrm{o}\mathrm{m}}+{\mathit{r}}_i\right)\times {\mathit{F}}_i\right),

where t, m_{\mathrm{b}}, {\mathit{I}}_{\mathrm{b}}, {\mathit{a}}_{\mathrm{b}}, {\text{ω}}_{\mathrm{b}}, and \mathit{g} are the time, body mass, body inertia, body linear acceleration vector, body angular velocity vector, and gravitational acceleration vector, respectively.

4 Landing control method

4.1 Landing method framework

The landing method controls the FLLWR to absorb the impact energy in the leg compression stroke to realize buffer landing. The maximum torque of IDU can be controlled to be less than the peak torque. Four legs are coordinated in accordance with the states of the legs.

As shown in Fig.8, the landing method framework consists of three layers: hardware input layer, landing method layer, and hardware output layer. The hardware input layer processes the actual joint torque {\mathit{\tau }}_{\mathrm{a}\mathrm{j}}, the actual joint angles {\mathit{\theta }}_{\mathrm{a}\mathrm{j}}, the actual roll angle {\theta }_{\mathrm{a}\mathrm{r}}, the actual pitch angle {\theta }_{\mathrm{a}\mathrm{p}}, and the corresponding velocities {\dot{\mathit{\theta }}}_{\mathrm{a}\mathrm{j}}, {\dot{\theta }}_{\mathrm{a}\mathrm{r}}, and {\dot{\theta }}_{\mathrm{a}\mathrm{p}} from the torque sensors, encoders, and IMU. In the landing method layer, if the amplitude of the tiptoe force for each leg in the vertical direction exceeds a threshold during the landing period, which signifies the occurrence of a touching-ground event, then the leg state will be changed to the grounded phase from the aerial phase. The trajectory planning module gives the reference trajectory of the joint angle {\mathit{\theta }}_{\mathrm{r}\mathrm{j}}, the reference roll angle {\theta }_{\mathrm{r}\mathrm{r}}, the reference pitch angle {\theta }_{\mathrm{r}\mathrm{p}}, and the corresponding velocities {\dot{\mathit{\theta }}}_{\mathrm{r}\mathrm{j}}, {\dot{\theta }}_{\mathrm{r}\mathrm{r}}, and {\dot{\theta }}_{\mathrm{r}\mathrm{p}}. Simultaneously, the active stiffness coefficient K and the damping coefficient B are given by the compliance planning. The actual tiptoe position can be obtained in real time by using the forward kinematics. The stabilizer is used to control the body roll angle, the body pitch angle, and the body height by generating the reference joint torque {\mathit{\tau }}_{\mathrm{r}\mathrm{j}}. The 12 joint impedance controllers generate the command joint torque {\mathit{\tau }}_{\mathrm{c}\mathrm{j}} of each IDU by using Eq. (18):

Fig.8 Landing method framework consisting of three layers: hardware inputs, landing method, and hardware outputs.

Full size|PPT slide

(18){\mathit{\tau }}_{\mathrm{c}\mathrm{j}}-{\mathit{\tau }}_{\mathrm{r}\mathrm{j}}=K\left({\mathit{\theta }}_{\mathrm{r}\mathrm{j}}-{\mathit{\theta }}_{\mathrm{a}\mathrm{j}}\right)+B\left({\dot{\mathit{\theta }}}_{\mathrm{r}\mathrm{j}}-{\dot{\mathit{\theta }}}_{\mathrm{a}\mathrm{j}}\right).

In the hardware output layer, the drivers transform {\mathit{\tau }}_{\mathrm{c}\mathrm{j}} into the corresponding current to control the motors.

4.2 State switching strategy

As illustrated in Fig.9, the robot state consists of the leg state and the body state. The leg state includes aerial phase and grounded phase, and the body state includes landing preparation phase, synergy adjustment phase, and full-supporting phase. The robot reaches above the landing site and uses the engines to adjusts the velocity and attitude to meet the landing condition. The lander starts to descend under the lunar gravity after turning off the engines. Concurrently, each leg switches to the aerial phase, and the body switches to the landing preparation phase. The legs touch the ground successively owing to the irregular terrain after a period of descent. Making a quick response to touch-ground event by the actual tiptoe force is difficult because its change is relatively small owing to the low coefficients of stiffness and damping in the landing prepare phase (Section 4.4). The virtual tiptoe force {\mathit{F}}_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}} generated by high virtual coefficients of stiffness (K_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}}) and damping (B_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}}) replaces the actual tiptoe force to judge the touchdown and to increase the response sensitivity to touching-ground event. {\mathit{F}}_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}} can be calculated as follows:

Fig.9 State switching method during buffer landing.

Full size|PPT slide

(19){\mathit{\tau }}_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}}=K_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}}\left({\mathit{\theta }}_{\mathrm{r}\mathrm{j}}-{\mathit{\theta }}_{\mathrm{a}\mathrm{j}}\right)+B_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}}\left({\dot{\mathit{\theta }}}_{\mathrm{r}\mathrm{j}}-{\dot{\mathit{\theta }}}_{\mathrm{a}\mathrm{j}}\right),

(20){\mathit{F}}_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}}={\left({\mathit{J}}_{\mathrm{v}}{(\mathit{q})}^{\mathrm{T}}\right)}^{-1}\left({\mathit{\tau }}_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}}-\mathit{H}\right),

where {\mathit{\tau }}_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}} is the virtual joint torque.

When the touching-ground event occurs at the first time, the body switches to the synergy adjustment phase in which some legs are grounded phase and the others are aerial phase. When all legs switch to the grounded phase, the robot state switches to the full-supporting phase.

4.3 Trajectory planning

The trajectory planning scheme for leg 4 is illustrated in Fig.10, and the other legs also adopt the same planning. In the landing preparation phase, the world coordinate system is assumed to coincide exactly with the body coordinate system. Therefore, the body position in the world coordinate system is zero, and the tiptoes keep initial positions: x_{\mathrm{t}\mathrm{i}\mathrm{p}0}, y_{\mathrm{t}\mathrm{i}\mathrm{p}0}, and {\mathit{\text{z}}}_{\mathrm{t}\mathrm{i}\mathrm{p}0}. An active compression distance of d is planned to reduce the impact because the impact energy is mainly distributed in the vertical direction. The initial tiptoe positions of all legs used in the experiments are listed in Tab.1. In the synergy adjustment phase, if the leg state changes to the grounded phase from the aerial phase, the net output torques of the thigh and shank IDUs will become 0 N·m due to the zero tiptoe force state (Section 4.4), resulting in a passive shrink under the action of the downward body. This design has two advantages. The lengths of the grounded legs will become shorter continuously, and the lengths of the aerial legs retain their initial value, which will contribute to the good adaption to landing terrain. The grounded legs will exert minimal influence on the body position and attitude because their forces are close to zero. At the beginning of the full-supporting phase, the world coordinate frame is updated to coincide with the current body coordinate frame, and the actual tiptoe positions, including x_{\mathrm{t}\mathrm{i}\mathrm{p}1}, y_{\mathrm{t}\mathrm{i}\mathrm{p}1}, and {\mathit{\text{z}}}_{\mathrm{t}\mathrm{i}\mathrm{p}1}, are recorded as a referenced tiptoe position in the whole supporting phases. For the body, an active body compression in the z direction from initial position {\mathit{\text{z}}}_{\mathrm{b}\mathrm{o}\mathrm{d}\mathrm{y}0} to target position {\mathit{\text{z}}}_{\mathrm{b}\mathrm{o}\mathrm{d}\mathrm{y}1} is planned to reduce the peak torque. The real-time body interpolation trajectory {\mathit{\text{z}}}_{\mathrm{i}\mathrm{n}\mathrm{t}\mathrm{e}\mathrm{r}\mathrm{p}} can be obtained from the following equation:

Tab.1 Initial tiptoe positions before buffer landing

Leg	xtip0/m	ytip0/m	ztip0/m
1	0.788	0.788	−0.750
2	−0.788	0.788	−0.750
3	−0.788	−0.788	−0.750
4	0.788	−0.788	−0.750

Fig.10 Trajectory planning of the body and tiptoe of leg 4 during buffer landing.

Full size|PPT slide

(21){\mathit{\text{z}}}_{\mathrm{i}\mathrm{n}\mathrm{t}\mathrm{e}\mathrm{r}\mathrm{p}}=c_0+c_1\left(t-t_0\right)+c_2{\left(t-t_0\right)}^2+c_3{\left(t-t_0\right)}^3+c_4{\left(t-t_0\right)}^4+c_5{\left(t-t_0\right)}^5+c_6{\left(t-t_0\right)}^6+c_7{\left(t-t_0\right)}^7,

(22)\{\begin{array}{l}[\begin{array}{cccc}{c}_0& c_1& c_2& c_3\end{array}]=[\begin{array}{cccc}{\mathit{\text{z}}}_{\mathrm{b}\mathrm{o}\mathrm{d}\mathrm{y}0}& v_0& a_0/2& j_0/6\end{array}],\\ c_4=\left(210d-T\left(\left(30a_0-15a_1\right)T+\left(4j_0+j_1\right)T^2+120v_0+90v_1\right)\right)/\left(6T^4\right),\\ c_5=\left(-168d+T\left(\left(20a_0-14a_1\right)T+\left(2j_0+j_1\right)T^2+90v_0+78v_1\right)\right)/\left(2T^5\right),\\ c_6=\left(420d-T\left(\left(20a_0-39a_1\right)T+\left(4j_0+3j_1\right)T^2+216v_0+204v_1\right)\right)/\left(6T^6\right),\\ c_7=\left(-120d+T\left(\left(12a_0-12a_1\right)T+\left(j_0+j_1\right)T^2+60v_0+60v_1\right)\right)/\left(6T^7\right),\\ [\begin{array}{cc}T& d\end{array}]=[\begin{array}{cc}{t}_1-t_0& {\mathit{\text{z}}}_{\mathrm{b}\mathrm{o}\mathrm{d}\mathrm{y}1}-{\mathit{\text{z}}}_{\mathrm{b}\mathrm{o}\mathrm{d}\mathrm{y}0}\end{array}],\end{array}

where v_0, v_1, a_0, a_1, j_0, and j_1 are the initial velocity, target velocity, initial acceleration, target acceleration, initial jerk, and target jerk of body in the z direction, respectively, t_0, t_1, t, and T are the initial time, end time, current time, and total duration time, respectively, and c_i (i=1,2,...,7) is the coefficient of interpolation trajectory. The equation can produce a trajectory of flexible impact by setting the velocity, acceleration, and jerk at the beginning and end time to zero. It avoids computing the inverse of 8 by 8 matrix in real time. The angle of each joint can be obtained by the inverse kinematics after completing the trajectory planning of the body and tiptoes during the whole landing period.

4.4 Compliance planning

The mechanical property of the leg is determined by the stiffness and damping coefficients, which can be adjusted actively through the compliance control. The compliance planning method for thigh and shank joints is shown in Fig.11, and the side joint keeps the high coefficients of stiffness (K_2) and damping (B_3) all the time. All thigh and shank joints keep the low stiffness coefficients K_0 and low damping coefficients B_0 to avoid a huge impact on each leg at the touching-ground moment in the landing preparation phase. At the beginning of synergy adjustment phase, the thigh and shank joints of grounded legs smoothly switch to zero tiptoe force state by the transition 1, and those of the aerial legs still retain the initial values K_0 and B_0. The damping coefficients increase to a median value B_2 by the transition 2 to avoid the torque overload because the angular velocity of active joint may reach its maximum at the start of full-supporting phase. Given that the angular velocity decreases after the energy dissipation in the transition 2, the damping coefficient can reach the maximum B_3 by the transition 3 to release the motor torque potential and realize the energy quick dissipation. The stiffness coefficient reaches the maximum K_2 by the transition 2 because the elastic force increases slowly from zero without an amplitude step. The switch equation of stiffness and damping coefficients is written as follows:

Fig.11 Compliance planning of stiffness and damping coefficients for thigh and shank joints.

Full size|PPT slide

(23)\begin{array}{rl}co{e}_{\mathrm{i}\mathrm{n}\mathrm{t}\mathrm{e}\mathrm{r}\mathrm{p}}=& -3c_0{\left(1-t_{\mathrm{r}}\right)}^2-6c_1\left(1-t_{\mathrm{r}}\right)t_{\mathrm{r}}+3c_1{\left(1-t_{\mathrm{r}}\right)}^2\\ & -3c_2{t}_{\mathrm{r}}^2+6c_2\left(1-t_{\mathrm{r}}\right)t_{\mathrm{r}}+3c_3{t}_{\mathrm{r}}^2,\end{array}

(24)\{\begin{array}{l}[\begin{array}{cc}{c}_0& c_3\end{array}]=[\begin{array}{cc}co{e}_0& co{e}_1\end{array}],\\ c_1=Tcoe{V}_0/3+co{e}_0,\\ c_2=-Tcoe{V}_1/3+co{e}_1,\\ T=t_1-t_0,\\ t_{\mathrm{r}}=\left(t-t_0\right)/T,\end{array}

where co{e}_{\mathrm{i}\mathrm{n}\mathrm{t}\mathrm{e}\mathrm{r}\mathrm{p}}, co{e}_0, co{e}_1, coe{V}_0, coe{V}_1, and t_{\mathrm{r}} are the current interpolation coefficient, initial coefficient, target coefficient, initial coefficient velocity, target coefficient velocity, and duration ratio time (from 0 to 1), respectively. The equation avoids solving the inverse of 4 by 4 matrix due to the explicit expression of c_0, c_1, c_2, and c_3.

4.5 Height and attitude stabilizer design

The resultant force in the vertical direction, including the leg force, the gravity force, and the inertial force, should be kept zero all the time in the full-supporting phase to maintain a stable height. When the roll and pitch angles of the body are not equal to the reference value, an extra virtual torsion is needed to control the body attitude. In accordance with the dynamic model, we have

(25)m_{\mathrm{b}}\mathit{g}+\sum _{i=1}^4{\mathit{F}}_{i\mathit{\text{z}}}=m_{\mathrm{b}}{\ddot{\mathit{z}}}_{\mathrm{b}},

(26)\sum _{i=1}^4\left(\left({\mathit{r}}_{\mathrm{c}\mathrm{o}\mathrm{m}}+{\mathit{r}}_i\right)\times {\mathit{F}}_{i\mathit{\text{z}}}\right)={\mathit{I}}_{\mathrm{b}}{\dot{\text{ω}}}_{\mathrm{b}},

where {\ddot{\mathit{z}}}_{\mathrm{b}} and {\dot{\text{ω}}}_{\mathrm{b}} are the body linear and angular acceleration, and {\mathit{F}}_{i\mathit{\text{z}}} is the ith leg force in vertical direction.

Here, the classical proportional–integral–derivative controllers are utilized to obtain the desired acceleration magnitude ({\ddot{\mathit{\text{z}}}}_{\mathrm{b}}) and angular acceleration {\dot{\text{ω}}}_{\mathrm{b}}=[\begin{array}{cc}{\dot{\omega }}_{\mathrm{b}x}& {\dot{\omega }}_{\mathrm{b}y}\end{array}{]}^{\mathrm{T}} as follows:

(27){\ddot{\mathit{\text{z}}}}_{\mathrm{b}}=k_{\mathrm{p},\mathit{\text{z}}}\left({\mathit{\text{z}}}_{\mathrm{r}\mathrm{b}}-{\mathit{\text{z}}}_{\mathrm{a}\mathrm{b}}\right)+k_{\mathrm{i},\mathit{\text{z}}}\int \left({\mathit{\text{z}}}_{\mathrm{r}\mathrm{b}}-{\mathit{\text{z}}}_{\mathrm{a}\mathrm{b}}\right)\mathrm{d}t+k_{\mathrm{d},\mathit{\text{z}}}\left({\dot{\mathit{\text{z}}}}_{\mathrm{r}\mathrm{b}}-{\dot{\mathit{\text{z}}}}_{\mathrm{a}\mathrm{b}}\right),

(28)\left[\begin{array}{c}{\dot{\omega }}_{\mathrm{b}x}\\ {\dot{\omega }}_{\mathrm{b}y}\end{array}\right]=k_{\mathrm{p},\theta }\left[\begin{array}{c}{\theta }_{\mathrm{r}\mathrm{r}}-{\theta }_{\mathrm{a}\mathrm{r}}\\ {\theta }_{\mathrm{r}\mathrm{p}}-{\theta }_{\mathrm{a}\mathrm{p}}\end{array}\right]+k_{\mathrm{i},\theta }\left[\begin{array}{c}\int \left({\theta }_{\mathrm{r}\mathrm{r}}-{\theta }_{\mathrm{a}\mathrm{r}}\right)\mathrm{d}t\\ \int \left({\theta }_{\mathrm{r}\mathrm{p}}-{\theta }_{\mathrm{a}\mathrm{p}}\right)\mathrm{d}t\end{array}\right]+k_{\mathrm{d},\theta }\left[\begin{array}{c}{\dot{\theta }}_{\mathrm{r}\mathrm{r}}-{\dot{\theta }}_{\mathrm{a}\mathrm{r}}\\ {\dot{\theta }}_{\mathrm{r}\mathrm{p}}-{\dot{\theta }}_{\mathrm{a}\mathrm{p}}\end{array}\right],

where k_{\mathrm{p},\mathit{\text{z}}}, k_{\mathrm{i},\mathit{\text{z}}}, and k_{\mathrm{d},\mathit{\text{z}}} are the proportional, integral, and derivative gain of the body in vertical direction, respectively, {\mathit{\text{z}}}_{\mathrm{r}\mathrm{b}}, {\mathit{\text{z}}}_{\mathrm{a}\mathrm{b}}, {\dot{\mathit{\text{z}}}}_{\mathrm{r}\mathrm{b}}, and {\dot{\mathit{\text{z}}}}_{\mathrm{a}\mathrm{b}} are the reference position, actual position, reference velocity, and actual velocity of the body in vertical direction, respectively, and k_{\mathrm{p},\theta }, k_{\mathrm{i},\theta }, and k_{\mathrm{d},\theta } are the proportional, integral, and derivative gains of roll and pitch angles.

Rewriting Eqs. (25)‒(28) in matrix form yields

(29)\begin{array}{rl}& \left[\begin{array}{cccc}1& 1& 1& 1\\ r_{\mathrm{c}\mathrm{o}\mathrm{m}y}+r_{1y}& r_{\mathrm{c}\mathrm{o}\mathrm{m}y}+r_{2y}& r_{\mathrm{c}\mathrm{o}\mathrm{m}y}+r_{3y}& r_{\mathrm{c}\mathrm{o}\mathrm{m}y}+r_{4y}\\ r_{\mathrm{c}\mathrm{o}\mathrm{m}x}+r_{1x}& r_{\mathrm{c}\mathrm{o}\mathrm{m}x}+r_{2x}& r_{\mathrm{c}\mathrm{o}\mathrm{m}x}+r_{3x}& r_{\mathrm{c}\mathrm{o}\mathrm{m}x}+r_{4x}\end{array}\right]\left[\begin{array}{c}{F}_{1\mathit{\text{z}}}\\ F_{2\mathit{\text{z}}}\\ F_{3\mathit{\text{z}}}\\ F_{4\mathit{\text{z}}}\end{array}\right]\\ =& \left[\begin{array}{c}{m}_{\mathrm{b}}\left({\ddot{\mathit{\text{z}}}}_{\mathrm{b}}-g\right)\\ I_{\mathrm{b}x}{\dot{\omega }}_{\mathrm{b}x}\\ I_{\mathrm{b}y}{\dot{\omega }}_{\mathrm{b}y}\end{array}\right].\end{array}

where g and F_{iz} are the third component of \mathit{g} and {\mathit{F}}_{iz}, respectively, r_{\mathrm{c}\mathrm{o}\mathrm{m}x}, r_{\mathrm{c}\mathrm{o}\mathrm{m}y}, r_{ix}, r_{iy}, r_{iz}, I_{\mathrm{b}x}, and I_{\mathrm{b}y} are the corresponding components of {\mathit{r}}_{\mathrm{c}\mathrm{o}\mathrm{m}}, {\mathit{r}}_i, and {\mathit{I}}_{\mathrm{b}}.

Equation (29) is an indefinite equation and has many feasible solutions. Researchers have proposed many effective methods to compute the Moore–Penrose inverse matrix, and some other researchers have obtained the tiptoe forces by solving a quadratic optimization problem [17,18]. Here, we propose a new method called virtual three-leg supporting algorithm (VTLSA) to obtain the solution of [\begin{array}{cccc}{F}_{1\mathit{\text{z}}}& F_{2\mathit{\text{z}}}& F_{3\mathit{\text{z}}}& F_{4\mathit{\text{z}}}\end{array}{]}^{\mathrm{T}}. As shown in Fig.12, the four-legged supporting problem is transformed into four groups of three-legged supporting model: leg 1-2-3 supporting, leg 1-2-4 supporting, leg 1-3-4 supporting, and leg 2-3-4 supporting. Each group is expected to provide a quarter of virtual resultant force in the z direction (F_{\mathrm{v}\mathrm{r}}) and virtual resultant torsion in the x and y directions ({\mathit{\tau }}_{\mathrm{v}\mathrm{r}}). When the virtual k\text{-}m\text{-}n leg supports the lander, we can obtain an analytical solution for the supporting force [\begin{array}{ccc}{F}_{k\mathit{\text{z}}}^{kmn}& F_{m\mathit{\text{z}}}^{kmn}& F_{n\mathit{\text{z}}}^{kmn}\end{array}{]}^{\mathrm{T}} from Eq. (30):

Fig.12 Virtual three-leg supporting algorithm.

Full size|PPT slide

(30)\left[\begin{array}{ccc}1& 1& 1\\ r_{\mathrm{c}\mathrm{o}\mathrm{m}y}+r_{ky}& r_{\mathrm{c}\mathrm{o}\mathrm{m}y}+r_{my}& r_{\mathrm{c}\mathrm{o}\mathrm{m}y}+r_{ny}\\ r_{\mathrm{c}\mathrm{o}\mathrm{m}x}+r_{kx}& r_{\mathrm{c}\mathrm{o}\mathrm{m}x}+r_{mx}& r_{\mathrm{c}\mathrm{o}\mathrm{m}x}+r_{nx}\end{array}\right]\left[\begin{array}{c}{F}_{k\mathit{\text{z}}}^{kmn}\\ F_{m\mathit{\text{z}}}^{kmn}\\ F_{n\mathit{\text{z}}}^{kmn}\end{array}\right]=\frac{1}{4}\left[\begin{array}{c}{m}_{\mathrm{b}}\left({\ddot{\mathit{\text{z}}}}_{\mathrm{b}}-g\right)\\ I_{\mathrm{b}x}{\dot{\omega }}_{\mathrm{b}x}\\ I_{\mathrm{b}y}{\dot{\omega }}_{\mathrm{b}y}\end{array}\right].

We can then obtain the expected leg force F_{i\mathit{\text{z}}} by summing the ith vertical force in leg k-m-n supporting (F_{i\mathit{\text{z}}}^{kmn} (i=k, m, or n)) for each group, which can be expressed as follows:

(31)\{\begin{array}{c}{F}_{1\mathit{\text{z}}}=F_{1\mathit{\text{z}}}^{123}+F_{1\mathit{\text{z}}}^{124}+F_{1\mathit{\text{z}}}^{134},\\ F_{2\mathit{\text{z}}}=F_{2\mathit{\text{z}}}^{123}+F_{2\mathit{\text{z}}}^{124}+F_{2\mathit{\text{z}}}^{234},\\ F_{3\mathit{\text{z}}}=F_{3\mathit{\text{z}}}^{123}+F_{3\mathit{\text{z}}}^{134}+F_{3\mathit{\text{z}}}^{234},\\ F_{4\mathit{\text{z}}}=F_{4\mathit{\text{z}}}^{124}+F_{4\mathit{\text{z}}}^{134}+F_{4\mathit{\text{z}}}^{234}.\end{array}

In accordance with Eq. (14), the real-time adjust tiptoe force {\mathit{F}}_{\mathrm{t}\mathrm{i}\mathrm{p}} can be transformed to the reference joint torque {\mathit{\tau }}_{\mathrm{r}\mathrm{j}}.

Currently, the Moore–Penrose inverse matrix can be computed quickly by using matrix decomposition techniques. The method provided by Katsikis and Pappas [33] can be utilized to solve large sparse matrix. The algorithm from Courrieu [34] was based on Cholesky factorization. Toutounian and Ataei [35] presented a CGS-MPi algorithm based on the conjugate Gram–Schmidt process, which is relatively robust for large sparse and rank deficient matrix. Katsikis et al. [36] developed a method based on QR factorization, which functions for singular square matrix or rectangular matrix. An efficient iterative algorithm provided by Esmaeili et al. [37] always converges with fourth-order. Tab.2 lists the computation time of six different algorithms and our method to obtain the tiptoe force with a tolerance of 10⁻¹³. For the Moore–Penrose inverse computation problem of small rectangular matrix that occurs during buffer landings, the methods of generalized inverse (Ginv) [33], tensor product matrix (TPM) [38], and improved Qrginv (IMqrg) [39], have the fastest solution speed. The methods of singular value decomposition (SVD), QR generalized inverse (Qrg) [36], and VTLSA have medium speed, and generalized inverse (Geninv) [34] method has the slowest solution speed. However, the difference among them is small, so they can all be used for stabilizer design of the FLLWR.

Tab.2 Computation time of different algorithms to solve Eq. (29)

Method	Time/s
SVD	3.81 × 10−5
Ginv	1.66 × 10−5
TPM	1.26 × 10−5
Geninv	5.10 × 10−5
Qrg	3.57 × 10−5
IMqrg	2.22 × 10−5
VTLSA	3.34 × 10−5

5 Landing numerical simulation

The numerical landing simulations based on the above kinematic and dynamic models are conducted by using the cosimulation of MATLAB and NX to verify the effectiveness of the landing control method and the safety of the driving system. The same parameters, such as the gravitational accelerations on the Moon g_{\mathrm{m}}=1.633 m/s², the ground stiffness coefficient K_g = 12000 N/m, and the ground damping coefficient B_g = 2000 N·s/m, are used in different simulations. The buffer landing results on the lunar surface are obviously different under different initial landing conditions: the vertical velocity v_{\mathit{\text{z}}}, the horizontal velocity v_{\mathrm{h}}, and the load mass. Theoretically, if the lander reaches above the landing site and the landing conditions are satisfied, the robot will finish the landing with a vertical impact. Hence, a vertical landing simulation is designed. When the lander descends from a height of 2.25 m with a total weight of 180 kg, that is, the load weight is 3.5 times of the self-weight, the vertical landing velocity reaches up to v_{\mathit{\text{z}}}=-2.1\mathrm{m}/\mathrm{s} at the touching-ground moment. A parabolic landing simulation (v_{\mathrm{h}}\ne 0 m/s) is designed to analyze the landing capacity against a horizontal disturbance velocity. Considering the symmetric distribution of the legs, v_x=v_{\mathrm{h}} is set to reduce the DoF of the lunar gravity testing platform from six to five. The initial conditions in the numerical simulations are listed in Tab.3. These parameters are selected as an example because they are typical conditions. However, other similar parameters can be used.

Tab.3 Initial conditions of numerical simulation

Situation	Weight/kg	vx/(m·s−1)	vz/(m·s−1)
Vertical landing	180	0.0	−2.1
Parabolic landing	180	0.4	−2.1

The results of vertical landing simulation are shown in Fig.13. At t=1.367\mathrm{s}, the torques of all legs change intensely, which signifies the occurrence of touching-ground, and the robot velocity in the z direction reaches up to the maximum −2.1 m/s. The legs then begin to compress and absorb the landing impact energy. The torque of the side IDU is about zero because the ground impact energy is mainly distributed in the vertical direction rather than the horizontal direction during the landing. At t=1.578\mathrm{s}, the peak torque reaches 155 N·m, and the torque margin is 22.5% relative to the rated peak torque 200 N·m. The thigh IDU torques T_{\mathrm{t}} and the shank IDU torques T_{\mathrm{s}} reach to the stable value after a damping oscillation period. At t=1.708\mathrm{s}, v_{\mathit{\text{z}}} is zero, and the robot body reaches the lowest height of 0.3983 m. During the whole landing process, the roll and pitch angles vary in a small range close to zero.

Fig.13 Results of vertical landing simulation: (a) joint torque in each leg, (b) body angle of roll and pitch, and (c) position and velocity of the body.

Full size|PPT slide

Fig.14 illustrates the results of the parabolic landing simulation. At t=1.367\mathrm{s}, all legs touch the ground, and the robot states are the same as the ones in the vertical landing simulation, except a horizontal disturbance velocity v_x=0.4 m/s. At t=1.635 s, the thigh IDU reaches up to peak torque 167.6 \mathrm{N}\cdot \mathrm{m} with a torque margin 16.2%. T_{\mathrm{t}} and T_{\mathrm{s}} reach the stable value after a damping oscillation period. At t=1.684\mathrm{s}, the robot body reaches the lowest height of 0.3965 m with v_z=0 m/s. Unlike the roll angle that remains in the initial value zero, the pitch angle fluctuates obviously with a maximum value of −7.898°, but it is quickly stabilized to be zero in about 3 s. Compared with the vertical landing simulation, the maximum torque is increased by 12.6 N·m, and the torque margin is decreased by 6.3% correspondingly.

Fig.14 Results of parabolic landing simulation: joint torques in (a) legs 1 and 4 and (b) legs 2 and 3; (c) body angle of roll and pitch; and (d) position and velocity of the body.

Full size|PPT slide

6 Prototype and experiments

6.1 Lunar gravity simulation platform

As shown in Fig.15, a 5-DoF-LGTP, including x, z, R_x, R_y, and R_{\mathit{\text{z}}} DoFs, is built to simulate the lunar gravity environment. The actual gravity of lander is the same as that on the Moon by using a counterweight system (Fig.16). This platform consists of two pillars, a horizontal guide rail, a cantilever, a wire rope, two fixed pulleys, two counterweights, and sandy soils. A laser tracker is used as the measurement equipment to record the position and velocity variations, and the angle variations are measured by an IMU. The mass of counterweights can be calculated by using the following equations:

Fig.15 Components of the 5-DoF-LGSP.

Full size|PPT slide

Fig.16 Counterweight system design of the 5-DoF-LGSP.

Full size|PPT slide

(32)\{\begin{array}{c}\left(m_{\mathrm{b}}+m_{\mathrm{c}2}\right)g_{\mathrm{e}}-m_{\mathrm{c}1}{g}_{\mathrm{e}}=\left(m_{\mathrm{b}}+m_{\mathrm{c}1}+m_{\mathrm{c}2}\right)g_{\mathrm{m}},\\ m_{\mathrm{b}}+m_{\mathrm{c}1}+m_{\mathrm{c}2}=m_{\mathrm{t}},\end{array}

where m_{\mathrm{b}}=40 kg is the body mass, m_{\mathrm{c}1} and m_{\mathrm{c}2} are the masses of counterweight 1 and counterweight 2, m_{\mathrm{t}} is the system mass of the lander and loads, and g_{\mathrm{m}} and g_{\mathrm{e}} are the gravitational accelerations on the Moon and Earth, respectively, that is, g_{\mathrm{m}}=g_{\mathrm{e}}/6. In this study, we aim to simulate the landing with a 140 kg payload, which indicates m_{\mathrm{t}}=140 kg. We obtain m_{\mathrm{c}1}=75\mathrm{k}\mathrm{g} and m_{\mathrm{c}2}= 65 kg.

6.2 Landing experiments

Similar to numerical simulations, two typical experiments, including the vertical landing and the parabolic landing, are conducted to validate the landing method and analyze the landing performance. The landing initial conditions are the same as those in the simulations.

As a representative example, Fig.17 shows the robot keyframe snapshots in the parabolic landing experiments. (1) is the initial position that is ready to land. (2) is the descending period in the air. (3) is the status at touching ground moment. (4) is compressing to the lowest position (5). (6) is extending to the highest position (7). The robot keeps the stable position after a damping oscillation (8).

Fig.17 Keyframe snapshots in the parabolic landing.

Full size|PPT slide

The vertical landing experiments are conducted, and the results are shown in Fig.18. The lander begins to land from a height of 2.25 m with a weight of 180 kg. At t=1.313\mathrm{s}, the torques have a sharp variation, and the lander touches the ground with a vertical velocity v_{\mathit{\text{z}}}=-2.1\mathrm{m}/\mathrm{s}. The lander vertical speed is zero, and the body is compressed to the lowest height of 0.3867 m after a deceleration for 0.192 s. The damping oscillation in the z direction lasts about 1 s, and then the body height is stabilized. The peak torque of the thigh IDU reaches up to −166.6 N·m with a margin of 16.7%. The pitch angle has a maximum deviation of −3.135°, but it is quickly adjusted to zero by the attitude stabilizer.

Fig.18 Results of vertical landing experiment: joint torques in (a) leg 1, (b) leg 2, (c) leg 3, and (d) leg 4; (e) body angle of roll and pitch; and (f) position and velocity of the body.

Full size|PPT slide

The parabolic landing experiments are conducted with the same condition as that in vertical landing experiments except for an extra horizontal disturbance velocity. The experiment results are illustrated in Fig.19. At t=1.46\mathrm{s}, the torque steps occur, signifying the ground-touching with a horizontal velocity v_x=0.4\mathrm{m}/\mathrm{s} and a vertical velocity v_{\mathit{\text{z}}}=2.1\mathrm{m}/\mathrm{s}. The vertical speed becomes zero, and the body is compressed to the lowest height of 0.3512 m after a deceleration for 0.189 s. The damping oscillation in the z direction lasts about 1.142 s, and then the body height is stabilized. The horizontal velocity is stabilized to zero after 1.893 s. The peak torque reaches up to −179.4 N·m with a margin of 10.3%. The roll angle is still about zero, and the pitch angle has a maximum deviation of −8.583°, but it is quickly controlled to be zero by the attitude stabilizer.

Fig.19 Results of parabolic landing experiment: joint torques in (a) leg 1, (b) leg 2, (c) leg 3, and (d) leg 4; (e) body angle of roll and pitch; and (f) position and velocity of the body.

Full size|PPT slide

The above experimental results are consistent with the simulation results. In the vertical landing, the small difference in the maximum torque between the simulation and the experiment is about 11.6 N·m, whereas it is about 11.8 N·m in the parabolic landing. The torque margin is sufficient with minimum value of 10.3% relative to the rated peak torque 200 N·m. In the simulations and experiments, the roll angle changes minimally with only 2° of fluctuation amplitude. The pitch angle changes significantly due to the horizontal impact from the x direction in the parabolic landing, reaching −7.898° and −8.583° in the simulation and experiment, respectively, but it is adjusted to zero in a short time through the height and attitude stabilizer. In the simulation and experiment, the lowest height parameter, which determines whether the engine nozzle at the bottom of the body collides with the ground, is higher than the design requirement (0.2575 m) that is determined by the distance from ground to the bottom of the body (0.2000 m) and the distance from the bottom of the body to the origin of the body (0.0575 m). Specifically, the lowest height in parabolic landing is lower than the one in vertical landing, and the minimum value 0.3512 m in parabolic landing still has a 36.4% safety margin. The proposed landing method is proved to be efficient by the successful landing in the simulations and ground experiments.

7 Conclusions and future work

In this study, an FLLWR is designed to execute the landing and roving missions on the Moon. To dissipate impact energy and achieve a safe landing, we propose a landing control method that consists of five parts: landing method framework, state switching strategy, trajectory planning, compliance planning, and height and attitude stabilizer. All sensor data, including encoder, force sensor, and IMU, are fused by the control method framework. The lander can make different responses to events in the landing process via the state switching. The height and attitude stabilizer are designed to realize the expected body height and attitude. A VTLSA is presented to quickly obtain the analytic solution of the stabilizer in the real-time control system. In the dynamic numerical simulation, the lander can achieve stable and safe landing under two typical working conditions: vertical landing and parabolic landing. The proposed control method is proved to be efficient for the landing in the lunar gravity by experiments on the 5-DoF-LGTP with a vertical landing velocity of ‒2.1 m/s and a payload of 140 kg, that is, 3.5 times of the self-weight. With or without an additional horizontal disturbance velocity of 0.4 m/s, the peak torque is 179.4 and 166.6 N·m, respectively, and the height safety margin, which avoids the collision between the engine nozzle installed on the bottom of the robot and the ground, is 36.4% and 50.1%, respectively.

The parameters satisfying the engineering application in compliance control are found by means of numerical dynamics simulation. In the future, we will study the relationship between landing and these parameters, and complete their optimal design systematically. The robot dynamic performance will be further explored to improve the buffer landing.

Nomenclature

Abbreviations
5-DoF-LGTP	Five-degree-of-freedom lunar gravity testing platform
COM	Center of mass
DoF	Degree-of-freedom
FLLWR	Four-legged landing and walking robot
Geninv	Generalized inverse based on Cholesky factorization
Ginv	Generalized inverse computation method
IDU	Integrated drive unit
IMqrg	Improved generalized inverse computation method based on QR factorization
IMU	Inertial measurement unit
LRV	Lunar roving vehicle
PFBM	Parallel five-bar mechanism
Qrg	Generalized inverse computation method based on QR factorization
RP	Revolute pair
SVD	Singular value decomposition
TPM	Tensor product matrix
VTLSA	Virtual three-leg supporting algorithm
Variables
a0, a1	Initial and target acceleration
ab	Body linear acceleration
B, B0, B1, B2, B3	Active damping coefficients
Bg	Ground damping coefficient
Bvirtual	Virtual damping coefficient
ci	Coefficient of interpolation trajectory
coeinterp, coe0, coe1, coeV0, coeV1	Current interpolation coefficient, initial coefficient, target coefficient, initial coefficient velocity, and target coefficient velocity in compliance planning, respectively
d	Active compression distance
Eps	Elastic potential energy of passive spring
Fi	ith tiptoe force from ground
Fiz, Fiz	Vector and the third component of the ith tiptoe force from ground in the z direction
F_{i\mathit{\text{z}}}^{kmn}	ith (i = k, m, or n) vertical force in leg k-m-n supporting
Ftip, Fx, Fy, Fz	Tiptoe force vector and its components
Fvirtual	Virtual tiptoe force
Fvr	Virtual resultant force in the z direction
g, g	Vector and the third component of gravitational acceleration
ge	Gravitational acceleration on the Earth
gm	Gravitational accelerations on the Moon
H	Correction vector
Ib, Ibx, Iby	Body inertia vector and its components
j0, j1	Initial and target jerk
Jv (q)	Velocity Jacobian matrix
K, K0, K1, K2	Active stiffness coefficients
kdz	Derivative gain in the z direction
kd,θ	Derivative gain of roll and pitch angles
Kg	Ground stiffness coefficient
ki,z	Integral gain in the z direction
ki,θ	Integral gain of roll and pitch angles
kp,z	Proportional gain in the z direction
kp,θ	Proportional gain of roll and pitch angles
Kps	Stiffness coefficient of passive spring
Kvirtual	Virtual stiffness coefficient
Ll, Lr	Distances from the intersection point of thigh and shank to the left and right fixed points of passive spring, respectively
Lps, Lps0	Current and original length of the passive spring, respectively
Ls	Shank length
Lt	Thigh length
mb	Body mass
mc1	Counterweight 1 mass
mc2	Counterweight 2 mass
mt	System mass of the lander and loads
Ob	Origin of body coordinate frame
Oli	Origin of the ith leg coordinate frame
Ow	Origin of world coordinate frame
{\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}}, {\dot{\mathit{P}}}_{\mathrm{t}\mathrm{i}\mathrm{p}}	Tiptoe position and velocity, respectively
{}^{\mathrm{b}}{\mathit{P}}_{\mathrm{l}i}	Origin position of {\mathrm{\Sigma }}_{\mathrm{l}i} in the body coordinate frame {\mathrm{\Sigma }}_{\mathrm{b}}
{}^{\mathrm{l}i}{\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}}	Tiptoe position in the leg coordinate frame {\mathrm{\Sigma }}_{\mathrm{l}i}
{}^{\mathrm{w}}{\mathit{P}}_{\mathrm{b}}	Body position in the world coordinate frame {\mathrm{\Sigma }}_{\mathrm{w}}
{}^{\mathrm{w}}{\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}}	Tiptoe position in the world coordinate frame {\mathrm{\Sigma }}_{\mathrm{w}}
\mathit{q}, \dot{\mathit{q}}	Generalized coordinate and velocity vector of joints, respectively
{\mathit{r}}_{\mathrm{c}\mathrm{o}\mathrm{m}}, r_{\mathrm{c}\mathrm{o}\mathrm{m}x}, r_{\mathrm{c}\mathrm{o}\mathrm{m}y}	Vector and its components from the COM to the origin Ob
{\mathit{r}}_i, r_{ix}, r_{iy}, r_{i\mathit{\text{z}}}	Vector and its components from the origin O_{\mathrm{b}} to the ith tiptoe
{}^{\mathrm{b}}{\mathit{R}}_{\mathrm{l}i}	Rotation transformation matrix from {\mathrm{\Sigma }}_{\mathrm{l}i} to {\mathrm{\Sigma }}_{\mathrm{b}}
{}^{\mathrm{w}}{\mathit{R}}_{\mathrm{b}}	Rotation transformation matrix from {\mathrm{\Sigma }}_{\mathrm{b}} to {\mathrm{\Sigma }}_{\mathrm{w}}
t, t_0, t_1, t_{\mathrm{r}}, T	Current time, initial time, end time, duration ratio time, and total duration time, respectively
T_{\mathrm{s}}, T_{\mathrm{t}}	Shank and thigh IDUs torques, respectively
v_0, v_1	Initial and target velocity, respectively
v_{\mathrm{h}}	Horizontal velocity
v_x	Horizontal velocity in the x direction
v_{\mathit{\text{z}}}	Vertical velocity
x_{\mathrm{b}}, y_{\mathrm{b}}, {\mathit{\text{z}}}_{\mathrm{b}}	Forward, left, and upper direction, respectively
x_{\mathrm{l}i}, y_{\mathrm{l}i}, {\mathit{\text{z}}}_{\mathrm{l}i}	Components of tiptoe position in the ith leg coordinate frame, respectively
x_{\mathrm{t}\mathrm{i}\mathrm{p}}, x_{\mathrm{t}\mathrm{i}\mathrm{p}0}, x_{\mathrm{t}\mathrm{i}\mathrm{p}1}	Current, initial position, and terminated position of tiptoe in the x direction, respectively
{\dot{x}}_{\mathrm{t}\mathrm{i}\mathrm{p}}, {\dot{y}}_{\mathrm{t}\mathrm{i}\mathrm{p}}, {\dot{\mathit{\text{z}}}}_{\mathrm{t}\mathrm{i}\mathrm{p}}	Current tiptoe velocity in the x, y and z directions, respectively
y_{\mathrm{t}\mathrm{i}\mathrm{p}}, y_{\mathrm{t}\mathrm{i}\mathrm{p}0}, y_{\mathrm{t}\mathrm{i}\mathrm{p}1}	Current, initial, and terminated positions of tiptoe in the y direction, respectively
{\mathit{\text{z}}}_{\mathrm{a}\mathrm{b}}, {\dot{\mathit{\text{z}}}}_{\mathrm{a}\mathrm{b}}	Actual body position and velocity in the z direction, respectively
{\ddot{\mathit{z}}}_{\mathrm{b}}, {\ddot{\mathit{\text{z}}}}_{\mathrm{b}}	Vector and the third component of body linear acceleration in the z direction, respectively
{\mathit{\text{z}}}_{\mathrm{b}\mathrm{o}\mathrm{d}\mathrm{y}0}, {\mathit{\text{z}}}_{\mathrm{b}\mathrm{o}\mathrm{d}\mathrm{y}1}	Body initial and target position in the z direction, respectively
{\mathit{\text{z}}}_{\mathrm{i}\mathrm{n}\mathrm{t}\mathrm{e}\mathrm{r}\mathrm{p}}	Real-time body interpolation trajectory
{\mathit{\text{z}}}_{\mathrm{r}\mathrm{b}}, {\dot{\mathit{\text{z}}}}_{\mathrm{r}\mathrm{b}}	Reference body position and velocity in the z direction, respectively
{\mathit{\text{z}}}_{\mathrm{t}\mathrm{i}\mathrm{p}}, {\mathit{\text{z}}}_{\mathrm{t}\mathrm{i}\mathrm{p}0}, {\mathit{\text{z}}}_{\mathrm{t}\mathrm{i}\mathrm{p}1}	Current, initial, and target positions of the tiptoe in the z direction, respectively
\text{δ}{\mathit{P}}_{\mathrm{t}\mathrm{i}\mathrm{p}}	Virtual displacement of tiptoe
\text{δ}\mathit{q}	Virtual displacement of the generalized coordinate vector \mathit{q}
\text{δ}{W}_{\mathrm{p}\mathrm{s}}	Virtual work of passive spring
{\mathit{\theta }}_{\mathrm{a}\mathrm{j}}, {\dot{\mathit{\theta }}}_{\mathrm{a}\mathrm{j}}	Actual joint angle and velocity, respectively
{\theta }_{\mathrm{a}\mathrm{p}}, {\dot{\theta }}_{\mathrm{a}\mathrm{p}}	Actual pitch angle and velocity, respectively
{\theta }_{\mathrm{a}\mathrm{r}}, {\dot{\theta }}_{\mathrm{a}\mathrm{r}}	Actual roll angle and velocity, respectively
{\theta }_{\mathrm{r}}, {\dot{\theta }}_{\mathrm{r}}	Rocker arm angle and velocity, respectively
{\mathit{\theta }}_{\mathrm{r}\mathrm{j}}, {\dot{\mathit{\theta }}}_{\mathrm{r}\mathrm{j}}	Reference joint angle and velocity, respectively
{\theta }_{\mathrm{r}\mathrm{p}}, {\dot{\theta }}_{\mathrm{r}\mathrm{p}}	Reference pitch angle and velocity, respectively
{\theta }_{\mathrm{r}\mathrm{r}}, {\dot{\theta }}_{\mathrm{r}\mathrm{r}}	Reference roll angle and velocity, respectively
{\theta }_{\mathrm{s}}, {\dot{\theta }}_{\mathrm{s}}	Side angle and velocity, respectively
{\theta }_{\mathrm{t}}, {\dot{\theta }}_{\mathrm{t}}	Thigh angle and velocity, respectively
\mathit{\tau }, {\tau }_{\mathrm{s}}, {\tau }_{\mathrm{t}}, {\tau }_{\mathrm{r}}	Joint torque vector and its components, respectively
{\mathit{\tau }}_{\mathrm{a}\mathrm{j}}, {\mathit{\tau }}_{\mathrm{c}\mathrm{j}}, {\mathit{\tau }}_{\mathrm{r}\mathrm{j}}	Actual, command, and reference joint torques, respectively
{\mathit{\tau }}_{\mathrm{v}\mathrm{i}\mathrm{r}\mathrm{t}\mathrm{u}\mathrm{a}\mathrm{l}}	Virtual joint torque
{\mathit{\tau }}_{\mathrm{v}\mathrm{r}}	Virtual resultant torsions in the x and y directions, respectively
{\mathrm{\Sigma }}_{\mathrm{b}}	Body coordinate frame
{\mathrm{\Sigma }}_{\mathrm{i}\mathrm{m}\mathrm{u}}	IMU coordinate frame
{\mathrm{\Sigma }}_{\mathrm{l}i}	ith leg coordinate frame
{\mathrm{\Sigma }}_{\mathrm{w}}	World coordinate frame
{\text{ω}}_{\mathrm{b}}	Body angular velocity
{\dot{\text{ω}}}_{\mathrm{b}}, {\dot{\omega }}_{\mathrm{b}x}, {\dot{\omega }}_{\mathrm{b}\mathrm{y}}	Vector and its components of body angular acceleration, respectively

Acknowledgements

This work was funded by the National Key R&D Program of China (Grant No. 2021YFF0307905).

Open Access

This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution, and reproduction in any medium or format as long as appropriate credit is given to the original author(s) and source, a link to the Creative Commons license is provided, and the changes made are indicated.

The images or other third-party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.

Visit http://creativecommons.org/licenses/by/4.0/ to view a copy of this license.

References

1 Lin R F , Guo W Z , Li M . Novel design of legged mobile landers with decoupled landing and walking functions containing a rhombus joint. Journal of Mechanisms and Robotics, 2018, 10(6): 061017

2 Lucas J W , Conel J E , Hagemeyer W A , Garipay R R , Saari J M . Lunar surface thermal characteristics from Surveyor 1. Journal of Geophysical Research, 1967, 72(2): 779–789

3 Williams R J , Gibson E K . The origin and stability of lunar goethite, hematite and magnetite. Earth and Planetary Science Letters, 1972, 17(1): 84–88

4 Gisler M , Sornette D . Exuberant innovations: the Apollo program. Society, 2009, 46(1): 55–68

5 Basilevsky A T , Abdrakhimov A M , Head J W , Pieters C M , Wu Y Z , Xiao L . Geologic characteristics of the Luna 17/Lunokhod 1 and Chang’e-3/Yutu landing sites, Northwest Mare Imbrium of the Moon. Planetary and Space Science, 2015, 117: 385–400

6 Ma Y Q , Liu S C , Sima B , Wen B , Peng S , Jia Y . A precise visual localisation method for the Chinese Chang’e-4 Yutu-2 rover. Photogrammetric Record, 2020, 35(169): 10–39

7 Xu Z , Guo D J , Liu J Z . Maria basalts chronology of the Chang’e-5 sampling site. Remote Sensing, 2021, 13(8): 1515

8 Zhang T , Zhang W M , Wang K , Gao S , Hou L , Ji J H , Ding X L . Drilling, sampling, and sample-handling system for China’s asteroid exploration mission. Acta Astronautica, 2017, 137: 192–204

9 Salzberg I M . Tracking the Apollo lunar rover with interferometry techniques. Proceedings of the IEEE, 1973, 61(9): 1233–1236

10 Kim D , Jorgensen S J , Lee J , Ahn J , Luo J W , Sentis L . Dynamic locomotion for passive-ankle biped robots and humanoids using whole-body locomotion control. The International Journal of Robotics Research, 2020, 39(8): 936–956

11 Raibert M , Blankespoor K , Nelson G , Playter R . BigDog, the rough-terrain quadruped robot. IFAC Proceedings Volumes, 2008, 41(2): 10822–10825

12 Semini C , Barasuol V , Goldsmith J , Frigerio M , Focchi M , Gao Y , Caldwell D G . Design of the hydraulically actuated, torque-controlled quadruped robot HyQ2Max. IEEE/ASME Transactions on Mechatronics, 2017, 22(2): 635–646

13 Bartsch S , Birnschein T , Cordes F , Kuehn D , Kampmann P , Hilljegerdes J , Planthaber S , Roemmermann M , Kirchner F . SpaceClimber: development of a six-legged climbing robot for space exploration. In: Proceedings of ISR 2010 (the 41st International Symposium on Robotics) and ROBOTIK 2010 (the 6th German Conference on Robotics). Munich: VDE, 2010, 1–8

14 Dirk S , Frank K . The bio-inspired scorpion robot: design, control & lessons learned. In: Zhang H X, ed. Climbing and Walking Robots: Towards New Applications. London: IntechOpen, 2007, 197–218

15 Bledt G , Powell M J , Katz B , Di Carlo J , Wensing P M , Kim S . MIT Cheetah 3: design and control of a robust, dynamic quadruped robot. In: Proceedings of 2018 IEEE/RSJ International Conference on Intelligent Robots and Systems. Madrid: IEEE, 2018, 2245–2252

16 Hutter M , Gehring C , Jud D , Lauber A , Bellicoso C D , Tsounis V , Hwangbo J , Bodie M , K R , Fankhauser S , P A , Bloesch M . ANYmal—a highly mobile and dynamic quadrupedal robot. In: Proceedings of 2016 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS). Daejeon: IEEE, 2016, 38–44

17 Arm P , Zenkl R , Barton P , Beglinger L , Dietsche A , Ferrazzini L , Hampp E , Hinder J , Huber C , Schaufelberger D , Schmitt F , Sun B , Stolz B , Kolvenbach H , Hutter M . SpaceBok: a dynamic legged robot for space exploration. In: Proceedings of 2019 International Conference on Robotics and Automation (ICRA). Montreal: IEEE, 2019, 6288–6294

18 Kolvenbach H , Hampp E , Barton P , Zenkl R , Hutter M . Towards jumping locomotion for quadruped robots on the moon. In: Proceedings of 2019 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS). Macao: IEEE, 2019, 5459–5466

19 Zhou J H , Chen M , Chen J B , Jia S . Optimal time-jerk trajectory planning for the landing and walking integration mechanism using adaptive genetic algorithm method. Review of Scientific Instruments, 2020, 91(4): 044501

20 Lin R F, Guo W Z, Li M, Hu Y, Han Y C. Novel design of a legged mobile lander for extraterrestrial planet exploration. International Journal of Advanced Robotic Systems, 2017, 14(6): 1729881417746120

21 Zhou J H , Jia S , Qian J C , Chen M , Chen J B . Improving the buffer energy absorption characteristics of movable lander-numerical and experimental studies. Materials, 2020, 13(15): 3340

22 Lin R F , Guo W Z , Zhao C J , Tang Y Y , Zhao C Y , Li Z Y . Topological design of a new family of legged mobile landers based on Truss-mechanism transformation method. Mechanism and Machine Theory, 2020, 149: 103787

23 Lin R F , Guo W Z . Type synthesis of reconfiguration parallel mechanisms transforming between trusses and mechanisms based on friction self-locking composite joints. Mechanism and Machine Theory, 2022, 168: 104597

24 Han Y C , Guo W Z , Peng Z K , He M D , Gao F , Yang J Z . Dimensional synthesis of the reconfigurable legged mobile lander with multi-mode and complex mechanism topology. Mechanism and Machine Theory, 2021, 155: 104097

25 Han Y C , Zhou C Z , Guo W Z . Singularity loci, bifurcated evolution routes, and configuration transitions of reconfigurable legged mobile lander from adjusting, landing, to roving. Journal of Mechanisms and Robotics, 2021, 13(4): 040903

26 Lorenz R D , Turtle E P , Barnes J W , Trainer M G , Adams D S , Hibbard K E , Sheldon C Z , Zacny K , Peplowski P N , Lawrence D J , Ravine M , McGee T G , Sotzen K , MacKenzie S M , Langelaan J , Schmitz S , Wolfarth L S , Bedini P D . Dragonfly: a rotorcraft lander concept for scientific exploration at Titan. Johns Hopkins APL Technical Digest, 2018, 34(3): 374–385

27 Balaram J , Aung M , Golombek M P . The Ingenuity helicopter on the perseverance rover. Space Science Reviews, 2021, 217(4): 56

28 Yin K , Zhou S L , Sun Q , Gao F . Lunar surface fault-tolerant soft-landing performance and experiment for a six-legged movable repetitive lander. Sensors, 2021, 21(17): 5680

29 Sun Q , Gao F , Chen X B . Towards dynamic alternating tripod trotting of a pony-sized hexapod robot for disaster rescuing based on multi-modal impedance control. Robotica, 2018, 36(7): 1048–1076

30 Luo J W , Wang S G , Zhao Y , Fu Y L . Variable stiffness control of series elastic actuated biped locomotion. Intelligent Service Robotics, 2018, 11(3): 225–235

31 Yin K , Sun Q , Gao F , Zhou S L . Lunar surface soft-landing analysis of a novel six-legged mobile lander with repetitive landing capacity. Proceedings of the Institution of Mechanical Engineers, Part C: Journal of Mechanical Engineering Science, 2022, 236(2): 1214–1233

32 McGhee R B , Frank A A . On the stability properties of quadruped creeping gaits. Mathematical Biosciences, 1968, 3: 331–351

33 Katsikis V N , Pappas D . Fast computing of the Moore–Penrose inverse matrix. Electronic Journal of Linear Algebra, 2008, 17: 637–650

34 Courrieu P . Fast computation of Moore–Penrose inverse matrices. Neural Information Processing—Letters and Reviews, 2005, 8(2): 25–29

35 Toutounian F , Ataei A . A new method for computing Moore–Penrose inverse matrices. Journal of Computational and Applied Mathematics, 2009, 228(1): 412–417

36 Katsikis V N , Pappas D , Petralias A . An improved method for the computation of the Moore–Penrose inverse matrix. Applied Mathematics and Computation, 2011, 217(23): 9828–9834

37 Esmaeili H , Erfanifar R , Rashidi M . An efficient method to compute the Moore–Penrose inverse. Advances in Pure and Applied Mathematics, 2018, 9(2): 143–152

38 Lu S X , Wang X Z , Zhang G Q , Zhou X . Effective algorithms of the Moore–Penrose inverse matrices for extreme learning machine. Intelligent Data Analysis, 2015, 19(4): 743–760

39 Ataei A . Improved Qrginv algorithm for computing Moore–Penrose inverse matrices. ISRN Applied Mathematics, 2014, 2014: 641706