Generalized Design, Modeling and Control Methodology for a Snake-like Aerial Robot

Snake-like robots have been developing in recent decades, and various bio-inspired ideas are deployed in both the mechanical and locomotion aspects. In recent years, several studies have proposed state-of-the-art snake-like aerial robots, which are beyond bio-inspiration. The achievement of snake-like aerial robots benefits both aerial maneuvering and manipulation, thereby having importance in various fields, such as industry surveillance and disaster rescue. In this work, we introduce our development of the modular aerial robot which can be considered a snake-like robot with high maneuverability in flight. To achieve such flight, we first proposed a unique thrust vectoring apparatus equipped with dual rotors to enable three-dimensional thrust force. Then, a generalized modeling method based on dynamics approximation is proposed to allocate the wrench in the center-of-gravity (CoG) frame to thrust forces and vectoring angles. We further developed a generalized control framework that can handle both under-actuated and fully actuated models. Finally, we show the experimental results with two different platforms to evaluate the flight stability of the proposed snake-like aerial robot. We believe that the proposed generalized methods can provide a solid foundation for the snake-like aerial robot and its applications regarding maneuvering and manipulation in midair.


Introduction
The development of snake-like robots has a long history that started with the study of real snakes' motion [1], which provided the mathematical basis. Then, various bioinspired mechanical designs were introduced to achieve the snake-like configurations [2][3][4][5], and the motion control and planning methods were then proposed for the snake-like locomotion [6][7][8][9][10]. The most unique feature of snake-like robots compared to other legged robots is obstacle-aided locomotion [11], which depends on explicit obstacles for pushing itself and moving forward. During the last decade, underwater snake-like robots have undergone considerable development [12][13][14]. These robots demonstrate a similar mechanism to obstacle-aided locomotion by utilizing the fin structure to generate sufficient reactive forces from the surrounding water.
The origins of terrestrial and underwater snake-like robots can be found in nature. However, it is significantly difficult to find phenomena related to snakes in the aerial domain because snakes have no wings, which are almost mandatory for flight. A unique case of snake motion in the air is via lateral undulation [15], where it uses the body as a wing to act like a glider. However, it is impossible to find a bio-inspired mechanism or motion to enable a snake-like robot to move freely in the air. Thus, we focused on the active propulsion mechanism (that is the propeller) to generate sufficient thrust force to handle the gravity. With this floating ability, the aerial robot can perform snake-like motions, such as squeezing through narrow spaces in midair, which can significantly extend the exploration range in various inspection and rescue applications. Furthermore, the ability of whole body grasping that mimics snake wrapping can be another advantage of a snake-like structure to handle and transport a large object in midair, indicating their superiority over ordinary aerial robots.
In our previous work [16], a unique propulsion system was presented to enable a snake-like aerial robot composed from four links. However, this modular design was not validated with different link numbers, especially with fewer links (e.g., two links) which showed an under-actuated property. Hence, in this work, we first developed a modularized link structure with the proper propulsion system for arbitrary link numbers, and further introduced generalized modeling and control methods to achieve stable flight for both under-actuated and fully actuated models, as shown in Figure 1.  Most snake-like robots have a highly modularized configuration, which contains at least one degree-of-freedom (DoF) at a joint to enable the articulated motion. The designs for snake-like robots can be categorized into six groups based on the type of locomotion: (i) robots with passive wheels [17]; (ii) robots with active wheels [18]; (iii) robots with active treads [4]; (iv) robots based on undulations using vertical waves [5]; (v) robots based on undulation using linear expansion [19]; (vi) robots with prismatic joints [20]. Most of the locomotion in these snake robots is obstacle-aided and thus non-holonomic, which indicates that the stability of the whole motion is easy to be achieved due to the large contact surface with the ground. Therefore, the active balance control, which is important for legged robots (particularly the bipedal robot), can be ignored for most of the snake-like robots.
The under water snake: the propulsion system, weak gravity environment and strong environment damping (water resistance). In the water domain, the locomotion depends on the balance between water resistance, buoyancy, and gravity. Thus, a fin structure is generally deployed to enable three-dimensional undulation locomotion in the water [21]. On the contrary, the screw as an active propulsion device is also used to achieve omnidirectional locomotion that makes it easier to perform obstacle avoidance and squeeze motion [22]. Given the sufficient damping effect from the surrounding water, the dynamics of underwater robots is relatively slow, which also indicates the easier controllability of the body balance compared to that of legged robots. However, the atmosphere has almost no resistance and buoyancy in most cases, which indicates that aerial robots are required to generate sufficient propulsion forces and a real-time control framework to keep their balance in the air.

Modular Aerial Robots
To achieve the snake-like motion in the air, a modular configuration is also necessary for the aerial robot. The aerial modular structure was first proposed in Distributed Flight Array [23], which introduces a modular flight unit composed of a single propeller. These flight units self-assemble on the ground before takeoff. Although arbitrary assembly shapes are available, the robot cannot change its shape in the air. Then, ModQuad was developed by [24] to enable self-assembly and self-disassembly in the air to provide the opportunity to change shape in the air. However, such a modular robot still has no ability to perform snake-like motion during a flight due to the lack of the joint structure.
Several modular aerial robots with joint connections were also proposed. For instance, a reconfigurable ModQuad connects modular quadrotors in a manner of a closed loop shape [25]. One of the effective applications of this closed-loop shape is to grasp objects using the inside ring structure. However, snake-like motion requires a chained configuration. Then, quadrotors with a chained configuration were proposed by [26], which shows the ability to squeeze narrow space using the snake-like motion. Similarly, a chained aerial robot that was composed of fully actuated flight units was developed by [27], which showed a high freedom of transformation in the air and also demonstrated aerial manipulation with the fixed-root link [28]. A common feature of these robots is that each flight unit has a complete flight ability, since there are more than four propellers in the unit. On the contrary, an articulated aerial robot proposed in [29] has a single propeller in each unit, which indicates the minimum actuators for snake-like motion in the air. However, the motion is only two-dimensional due to the under-actuated property.
To achieve greater maneuvering with joint motion in the air, a unique thrust vectoring apparatus with a dual-rotor was proposed in [16], which resulted in a snake-like aerial robot called DRAGON. This robot has two DoF in each joint module, which enables the highest freedom in both aerial maneuvering and manipulation [30,31]. However, this proposed mechanical design was only validated in four-link type, and the feasibility with fewer links that is an under-actuated model is unknown. Although a similar configuration was applied in the two-link model presented in [32], this robot only showed the stability with the fixed root link, and the flight stability has not been achieved yet. Therefore, the aim of this work was to obtain a generalized design, and modeling and control methods for arbitrary link numbers that can be either under-or fully actuated.

Contributions
Based on our previous works in [16,[30][31][32], we propose the generalized methodologies for both under-actuated and fully actuated models that can be summarized as follows: • We introduced a generalized design for snake-like aerial robots, including a thrust vectoring apparatus with dual rotors that can generate different thrust forces. • We presented a generalized modeling method for articulated aerial robots based on an approximated model and further proposed two different actuator allocation strategies according to the number of vectoring apparatus. • We developed a generalized control framework that utilizes the proposed actuator allocation to enable the stable flight for both under-actuated and fully actuated models. • We performed experiments with two different platforms, as shown in Figure 1, to demonstrate the feasibility of the proposed design, and modeling and control methods for the flight with joint motion in midair.
We did these things to validate our proposed modeling and control methods, and achieve stable flight for both two and four-link types.

Notation
All the symbols in this paper are explained upon first appearance. Boldface symbols (most are lowercase, e.g., r) denote vectors, whereas non-boldface symbols (e.g., m or I) denote either scalars or matrices. A coordinate regarding a vector or a matrix is denoted by a left superscript, e.g., {A} r expresses r with reference to (w.r.t.) the frame {A}. We define {W} as a unique frame to represent the inertial reference frame. Then, subscripts are used to express: a target frame for a vector or matrix; e.g., {W} r {A} represents the 3D position of the frame {A} w.r.t {W}-and/or a relation or attribute for a scalar; e.g., {W} r {A} x represents the scalar position of the frame {A} along the x axis of the frame {W}.

Organization
The remainder of this paper is organized as follows. The generalized design for snake-like aerial robot is described in Section 2. Then, the generalized modeling method is presented in Section 3, which is followed by the generalized control framework in Section 4. Finally, we show the flight experiments using an under-actuated and a fully actuated robot in Section 5. The conclusions are summarized in Section 6.

Generalized Design
In this section, we propose a generalized design of an articulated aerial robot with distributed rotors, as shown in Figure 2A. Here, the main link configuration of the proposed robot is cylindrical, which is connected by an actuated joint to enable the snake-like motion. The link does not necessarily contain the thrust vectoring apparatus (i.e., link3, link5, and link 7). (B) Two-DoF thrust vectoring apparatuses (θ j , φ j ). {G j } is a frame attached to the origin of vectoring apparatus, and the x axis is aligned with the x axis of {L i } and rotates around it with φ j . {F j1 } and {F j2 } are the frames attached to the dual rotors, and {F j } is a frame in the middle of the them, where the z axis is parallel to the rotor rotation axis and is titled from the z axis of {G j } with θ j . λ j1 and λ j2 are the thrust forces generated by the dual rotors. (C) Two-DoF joint module composed of two orthogonal joint axes (q i_yaw , q i_pitch ).

Dual-Rotor Vectoring Apparatus
For a stable flight of this robot, the design of the thrust apparatus is significantly important. We put effort on the design of the minimal thrust apparatus that enables stable flight for an arbitrary number of links, which is based on the following two aspects: the vectoring DoF and the rotor number.
For the vectoring DoF, we consider a 2-DoF vectoring apparatus to achieve the three dimensional thrust force. As shown in Figure 2B, this concept is achieved by equipping two perpendicular vectoring angles (φ j , θ j ). This allows the thrust force to point any direction in arbitrary posture and maximize its performance in the flight control. This pair of vectoring axes are controlled by two independent servos with a general PID position control to track the desired angles. On the other hand, θ i = ±90 • would cause a singularity, because the thrust force can only point along the direction of link rod, regardless of the change in φ i .
Next, we determine the number of rotors on a vectoring apparatus. To achieve the minimal configuration, the number of rotors should be as low as possible. For a single rotor deployment, the optimal rotor position should be {F j }, as shown in Figure 2B, which, however, would collide with the vectoring shaft. Furthermore, only a single rotor in each vectoring apparatus (that is in each link module) would also induce the singular configuration that cannot fly, such as the straight-line configuration. Therefore, we use a dual-rotor structure, as shown in Figure 2B, and allow the rotors to generate different thrusts. Then, the articulated robots with distributed rotors can fly in a straight-line configuration with more than two vectoring apparatuses. It is notable that the number of the vectoring apparatuses is not necessarily deployed in all links, and the relative position in each link is also arbitrary. Nevertheless, deployment at the center of the link can provide a better weight balance than other position.

Two-DoF Joint Module
To achieve a snake-like motion, the DoF of the link pose w.r.t. the neighboring link should be more than two, which implies two orthogonal joints are necessary between two links. Therefore, we introduce a composite 2-DoF joint module shown in Figure 2C that consists of two identical single joint structures and an inter connection part. The first joint rotates about the pitch axis that corresponds to the angle of q i_pitch , and the second joint rotates about the orthogonal yaw axis that corresponds to the angle of q i_yaw . Each single joint is actuated by a servo independently. Similarly to the servo for the vectoring apparatus, the general PID position control is used for joint motion.
Two typicals model based on the proposed design are depicted in Figure 3. From the next section on, we will discuss the dynamics and control for those models.

Generalized Modeling Method
In our work, we address the articulated model with more than two links. To achieve the aerial transformation by such a multilinked model, a comprehensive investigation on modeling is important. In this section, we first describe the approximation method to obtain the simplified multilinked model, and then present the actuator allocation for both under-actuated and fully actuated models.

Approximation Model
As shown in Figure 2, the kinematic model of the proposed aerial robot is composed from a chained link structure. We assume the number of links is N L ; then the vectoring of joint angles q ∈ R 2(N L −1) can be defined as q := q 1_yaw , q 2_pitch , q 2_yaw , q 2_pitch , · · · . The thrust vectoring apparatus, as shown in Figure 2B, consists of two vectoring angles (θ j , φ j ) and dual rotors that generate two thrust forces (λ j1 , λ j2 ). Therefore, there are four control inputs, θ j , φ j , λ j1 , and λ j2 , in each vectoring apparatus, and we developed two different usage patterns for the under-actuated and fully actuated models, respectively. It is also notable that the vectoring apparatus is not necessarily deployed in each link module, as shown in Figure 2A. Therefore, the number of rotors, N r , can be different from 2N L .
Then, the dynamic model of such multilinked model w.r.t. the entire CoG frame {CoG} can be written as follows: where the first equation denotes the dynamic motion of the entire linear momentum, which is described in the inertial frame {W}, whereas the second equation denotes the dynamic motion of the entire rotational momentum, which is described in the CoG frame of the entire multibody model (i.e., {CoG}). The third equation corresponds to the joint motion. g is a three-dimensional vector expressing gravity. {W} P ∑ and {CoG} L ∑ on the left sides of (1) and (2) are the total linear and angular momentum, respectively, which are both affected by the joint angles, vectoring angles, and their velocities, whereas {CoG} f and {CoG} τ on the right sides are the total wrench obtained from all vectored thrust forces. The allocation from the vectored thrust forces from this wrench is the key to achieving the flight control, which is described in Section 3.2 in detail.
In Equation (3), M J (q) denotes the inertial matrix, whereas c(q,q) is the term related to the centrifugal and Coriolis forces in joint motion. The symbol "s" stands for "segment" in multilinks. J r i ∈ R 3×N J and J s i ∈ R 3×N J are the Jacobian matrices for the frames of the i-th rotor and the i-th segment's CoG, respectively. τ q ∈ R N J is the vector of joint torque, and f i denotes the three dimensional force generated by each vectored rotor.
The entire dynamics model summarized in (1)∼(3) shows the high complexity due to the joint motion, and thus the real-time feedback control based on such a nonlinear model is significantly difficult. Therefore, a crucial quasi-static assumption is introduced in our work to simplify the dynamics; i.e., all the joints are actuated well and slowly by servos (q ≈ 0;q ≈ 0). Then, the joint velocity and acceleration can be omitted regardless of the joint motion. Under this assumption, the original dynamic model can be approximated as follows: Equations (4) and (5) still show the properties of the time-variant model because q changes over time and affects both the cog position {W} r {CoG} (q) and the overall rotational inertia I Σ (q). (6) shows the equilibrium between the joint torque, thrust force, and gravity, which can help us to obtain the desired joint torque from the thrust force. Given that we assume that joints are well controlled by the feedback position control of servos, it is indicated that there is no necessity to perform torque control based on (6). By ignoring (6), the whole model for control can be finally considered as a single rigid body by our approximation.

Actuator Allocation
Allocation from the three-dimensional thrust forces f i to the CoG wrench {CoG} w(:= [ {CoG} f {CoG} τ] T ) provides the connection from the whole body feedback control to the actuators that includes the thrust force λ and vectoring angles θ and φ.
The force gap (λ j1 − λ j2 ) from the dual rotors, as shown in Figure 2B, would induce a moment load on the servo that controls the vectoring angle φ. Since this vectoring servo should be compact and thus weak, it is considered difficult to dynamically control both the vectoring angle φ j and the thrust forces λ j1 and λ j2 at the same time.
Therefore, we developed two different strategies as follows: (1) a dual-rotor mode that allows different forces for λ j1 and λ j2 , but the vectoring angle φ j is constant in control framework; (2) a virtual-single-rotor mode that assigns the same force for λ j1 and λ j2 to avoid the moment load on vectoring angle φ i , and thus, we can use φ i as the control input.

Dual-Rotor Mode
The force and torque related to the j-th rotor module can be written as: where * denotes the index of the rotor, which is either 1 or 2. u j is the thrust unit normal, κ is the ratio of rotor thrust to its drag, and σ j * is the rotational direction of each rotor. {CoG} p {Fj * } (q, φ j ) is the rotor position that is affected by the joint angles q and the vectoring angle φ j . Then, the relationship between the target wrench and rotor thrust can be expressed by whereÑ r = N r 2 is the number of rotor apparatus, and {CoG} u j * , {CoG} v j * and λ j * correspond to (7) and (8).

Virtual-Single-Rotor Mode
Given that the dual rotors generate the same thrust forces, there is no moment that occurs in the vectoring angle φ. Then, we can count the pair of rotors as an integrated rotor that generates a combined uni-directional thrust λ j = λ j1 + λ j2 . In addition, the drag moment and gyroscopic moment can be ideally counteracted. Then, the force {CoG} f j and torque {CoG} τ j related to the j-th rotor module can be written as: where b 3 = [0 0 1] T . Definitions of the frames {L j }, {G j }, and {F j } can be found in Figure 2B. {CoG} p {Fj} in (12) is the position of the frame {F j }, which depends on the joint angles q and the vectoring roll angle φ j because there is an offset from {G j } to {F j }, as shown in Figure 2C.
Then, the total wrench in the CoG frame can be given by where {CoG} u j , {CoG} v j and λ j correspond to (11) and (12). If allocation matrix Q in (13) is full-rank, an arbitrary wrench {CoG} w can be achieved by the control input of λ j , φ j , and θ j , which can be considered fully actuated. If a model has more than two rotor vectoring apparatuses, the full pose control can be achieved by using this allocation mode. However, only two apparatuses imply the bijection between six control input (i.e., λ 1 , λ 2 , φ 1 , φ 2 , θ 1 , and θ 2 ) and full pose motion SE (3), which can easily result in a control input that exceeds the valid range, especially for the thrust force λ. However, for a model with more than three vectoring apparatuses, there is the redundancy in control input. Therefore, we apply the dual-rotor mode for model with two vectoring apparatuses and the virtual-single-rotor mode for other cases.

Common Framework
Based on the approximated model and the allocation strategy proposed in Section 3, we present a common framework for both under-actuated and fully actuated models, as shown in Figure 4, which contains the first part for the pose control in the entire CoG motion, which is followed by the control allocation with the allocation strategies proposed in Section 3.2. Figure 4. Overview of the control framework developed in this work, which is a part of the whole system. "Model approximation" is presented in Section 3.

Full Pose Control
For the approximated dynamics (4) and (5), feedback control based on a common PID control is applied as follows: where K f , * are the PID gain diagonal matrices.
The attitude control follows the SO(3) control method proposed by [33]: where [ ] ∨ is the inverse of a skew map, and R := {W} R {CoG} , ω := {CoG} ω for convenience. The efficiency of this full pose control ( (15) and (17)) for the fully actuated model has been validated in [31], and thus we extend the usage to the under-actuated model with a proper control allocation in this work.

Under-Actuated Model
To stabilize the attitude and altitude for the under-actuated model, a truncated desired wrench can be given byw A truncated allocation matrixQ from (21) is also introduced as follows: where {CoG} u z j * is the third element of {CoG} u j * . Then, the desired thrust forces λ des can be given by λ =Q #wdes (22) where (·) # denotes the weighted MP-pseudo-inverse. For under-actuated position control in the horizontal directions, the attitude in the roll and pitch directions is generally utilized. However, the proposed under-actuated robot can use both attitude and vectoring apparatus angles for x and y control. Note that we assume that these angles are sufficiently small, and these angles do not affect the attitude and altitude control. For the position control, the attitude in the roll and pitch angles can be expressed as follows: where (φ d , θ d , ψ d ) are the XYZ-Euler angles of the CoG frame orientation. For the proposed rotor-distributed robots, it is difficult to stabilize the flight in some configurations due to the large moment of inertia. To address this problem, we also use the i-th thrust vectoring angles θ des j , φ des j as follows: where n j is expressed by

Fully Actuated Model
The complete desired wrench can be summarized as follows: The control objective is to calculate the desired thrust λ des and the desired vectoring angles φ des , θ des from the desired CoG wrench w des from (29). For a model with more than three rotor vectoring apparatuses (i.e.,Ñ r > 2), there is infinite solution of (λ, φ, θ) according to (13). Then, an optimal geometry allocation can be given by where Q is derived from (14). Given that (31) is a nonlinear constraint for this optimization problem, the computational cost is relatively large, and the calculation time would increase asÑ r becomes larger. Therefore, we developed an iterative solution using a gradient to guarantee the convergence to the optimal (or at least suboptimal) solution. We first follow the method proposed by [34], which utilizes the vectored forces f j ∈ R 3 from (11) . Then, the above optimization problem can be rewritten as: where E 3×3 is a 3 × 3 identity matrix and [·×] denotes the skew symmetric matrix of a three-dimensional vector. Then, the closed-form for (32) and (33) and the desired thrust and vectoring angles can be directly given by: ).
Given the offset from the frame {G j } to frame {F j } in the two-DoF vectoring apparatus, as shown in Figure 2B, the result of vectoring angles φ and θ from (38) and (39) will change the position {CoG} p {Fi} in (35) again. Thus, the results of (37)∼(39) will no longer satisfy the constraint (33) becauseQ has changed again. Therefore, we introduce the residual term for (31) using the results of (37)∼(39), and also compute the derivative with respect to λ, φ, and θ: The partial derivative elements in (44) and (45) can be calculated from the multilinked kinematics model.
Our goal is to find a solution of λ, φ, and θ to guarantee zero . Then, we start from the initial values λ 0 , φ 0 , and θ 0 calculated from (37)∼(39), and perform the following linear iteration with the objective of minimizing ∑Ñ r j=1 δλ j 2 + δφ j 2 + δθ j 2 at each iteration: where J # w is the psuedo-inverse matrix of J w , and k ∈ [0, 1, 2, · · · ] is the iteration number. The most computationally intensive operation in this iteration process is the calculation of the inverse matrix with a size of the 6 × 6 for J # w , which can be solved instantaneously for real-time control. In most of the cases, it only requires 2 or 3 iterations to get a sufficiently small value of (i.e., ∼ 10 −6 ). This may be attributed to the initial values λ 0 , φ 0 , and θ 0 being relatively close to the convergent solution. We finally define the convergent solution as λ des , φ des , and θ des .

Experiments
In this section, we present the development of platforms based on the generalized design, modeling, and control methods, and further demonstrate the experimental results regarding the aerial transformation motion.

Platforms
We developed two different types of platforms: (1) an under-actuated robot with three links and two vectoring apparatuses and (2) a fully actuated robot with four links and four vectoring apparatuses. The main pipes of these robots were made of carbonfiber-reinforced plastic (CFRP) whose diameter was 25 mm and thickness was 1mm, and the power cables were in these pipes. Thus, these links were connected by a joint unit made of aluminum plates. The onboard system diagram is summarized in Figure 5. The robot had an onboard computer and a main control board to perform the state estimation, flight controller, and motion planning. It is notable that the process of flight control can be divided into two parts: (1) process with heavy computation (e.g., the matrix inversion (22) and the iteration process (46)∼(47)) that performed in the onboard computer; (2) real-time process (e.g., attitude control (17)) that performed in the main control board. Here, the attitude was estimated using IMU on this board at the rate of 500 Hz, and the estimated error was less than 1 • . In addition, an external motion capture system, of which the position estimated error was less than 1 mm, was applied in our experiment to obtain the position of the root link. Furthermore, each link was equipped with a small control board called neuron that connected with the main control board via control area network (CAN). The joints connecting these links were actuated by servo motors, of which the measurement error of joint angle was less than 0.1 • .

Motion Planner
Target Trajectory Figure 5. The onboard system diagram of the proposed snake-like aerial robots. The system is composed of three parts: (1) an onboard computer to perform the processes that require extensive computational resources; (2) main control board to perform real-time processes, such as the attitude estimate and control; (3) a neuron in each link module to transmit actuator commands from the main control board.

Under-Actuated Model
The detailed configuration of the under-actuated platform is depicted in Figure 6. The under-actuated robot was composed of three links (root link, middle link, and end link), and root and middle links were equipped with vectoring apparatuses. The main specifications of this robot can be found in Table 1. The CPU of the onboard PC was an Intel m3-8100Y with quad cores. The joint servo motors were XH430-W350-R (Dynamixel, stall torque: 3.4 Nm), and the vectoring apparatus joints for φ and θ axes used XH430-W350-R and XL430-W250-T((Dynamixel, stall torque: 1.5 Nm)), respectively. Here, two vectoring apparatuses were the minimal configuration for the articulated robot flight. Furthermore, the diameter of the three-blade propeller was 5 inch, and the maximum thrust was up to 20 N when the voltage was 25.2 V.

Vectoring apparatus Onboard Computer
Link joint

Root vectoring apparatus
Root Link

Middle Link
End Link Yaw joint Pitch joint Figure 6. Under-actuated robot platform with distributed rotors. This robot was composed of three links and two vectoring apparatuses, and the root and middle links were equipped with the vectoring apparatus. The two rotors on the vectoring apparatus generated different levels of thrust.

. Fully Actuated Model
A fully actuated platform was composed of four link modules with four rotor vectoring apparatus as shown in Figure 7. The main specifications of this robot can be found in Table 2. Furthermore, processes related to modeling, control, and motion planning were all performed inside the on-board compact computer, of which the CPU was an Intel Atom x7-Z8700 with quad cores. The main difference from the prototype developed in our previous work [16] was the propulsion system (i.e., rotor and propeller). We improved the thrust performance by increasing the blades number from 6 to 12. More blades resulted in slower rotation speed. Furthermore, we chose an inner rotor that has a KV rate of 2100 KV, which can significantly suppress the vibration due to the rotor high-speed rotation (i.e., 20,000∼30,000 RPM). Onboard batteries were deployed for the flight. Two LiHv batteries (1300 mAh, 22.8 V) were attached at each link to provide power for a pair of rotors, which enabled the maximum flight time of 3 min.

Pitch joint
Yaw joint Onboard Computer

Flight Experiments
The main goal of the experiments in this work is to verify the feasibility of the proposed generalized modeling and control methods for both under-actuated and fully actuated models that are the crucial foundation to performing complex applications. Therefore, we focus on evaluation of the flight stability and robustness with a given joint motion.

Under-Actuated Robot
Aerial maneuvering with joint motion was evaluated with the under-actuated threelink robot. The configuration was changed during flight, as shown in Figure 8A. The under-actuated flight control framework, as described in Section 4.3.1, was used, and the control parameters used in the experiments are summarized in Table 3. The position control ran at 40 Hz on the onboard processor, whereas the attitude control ran at 200 Hz on the main control board.
0.0, 0.0 Figure 8B-E plot the trajectories of thrust forces, joint angles, position, and attitude errors, respectively. In this experiment, the trajectories of thrust forces did not change significantly during hovering with joint motion, as shown in Figure 8B Table 4. During the configuration change, the flight stability could still be guaranteed, which demonstrated the feasibility of the proposed modeling and control methods for under-actuated model. To evaluate the flight performance of the proposed snake-like aerial robot in high places, we further conducted an experiment that involved large elevation, as shown in Figure 9A. During this experiment, the robot ascended from the floor to the ceiling at the height of approximately 3.0 m, as plotted in Figure 9B. The position-tracking errors during the whole flight are plotted in Figure 9C; the RMSE were [0.099, 0.035, 0.058] m. This result not only verified the stability of the proposed snake-like aerial during the ascending motion, but also demonstrated an efficiency in the elevated terrain that is close to the ceiling.

Fully Actuated Robot
A large-scale aerial maneuvering with joint motion was evaluated with the fully actuated platform. Based on the motion planning method proposed in [30], we designed a joint trajectory that can squeeze a small opening in midair like a snake, as depicted in Figure 10. However, in the actual experiment, we omitted the opening and ceiling with the aim of providing a better visualization for the trajectory tracking, as shown in Figure 11. The fully actuated flight control framework, as described in Section 4.3.2, was used, and the control parameters used in the experiments are summarized in Table 5. Both the motion planning and flight control frameworks were performed on the onboard processor.   Figure 10. Joint motion to squeeze a small opening in midair. The width of opening is smaller than that of the robot under the normal form. In addition, a low ceiling also limits the free space. Therefore, a snake-like squeezing motion is required. Figure 11. Complex maneuvering by the proposed snake-like fully actuated aerial robot in midair.  plot the trajectories of the CoG pose, joint angles, rotor vectoring angles, and thrust forces, respectively. The trajectory of the z-axis of the CoG motion indicates an ascending motion, whereas other axes show more complex motion. Furthermore, the joints changed in a large range (i.e, [−1.0, 1.6]), as shown in Figure 13, indicating the feasibility of complex maneuvering by our proposed robot platform. The maximum joint velocity was 0.31 rad/s. Although this velocity broke the quasi-static assumption that we applied in the modeling method, the flight stability during the joint motion was still guaranteed, which showed the robustness of our flight method for the joint motion in midair. The trajectories of the rotor vectoring angles (φ, θ) in Figure 15 and the thrust forces λ in Figure 15 showed the behavior of the outputs from the proposed flight framework, which contains the fluctuations that resulted from the D control term in both (15) and (16).  The RMSE of position and rotation control during the whole motion is summarized in Table 6. These results demonstrate the relatively high accuracy of full-pose tracking in the situation that involved joint motion, which indicated the feasibility of the modeling and control methods proposed in Sections 3 and 4. Regarding the robustness against the external disturbance, we considered the ground effect. During the joint motion, the minimum distance from the floor to the lowest link was less than 0.2 m, which corresponds to the time of 4 in Figure 11. The lowest link had significantly higher interference caused by the downwash compared to other higher links, and this interference can be considered as an external force acting at the lowest link. Nevertheless, Figure 12 still demonstrated a relatively small deviation from the desired trajectory of the CoG pose, which indicated that our proposed control method can guarantee promising robustness against varying external disturbance.   The high trackability for the complex trajectory, as shown in Figures 12 and 13, also indicated the potential to perform snake-like squeezing motion in midair, which can benefit the exploration in tight and elevated terrains, such as inspection between the plant pipelines in high place.
In terms of design, the above experiments demonstrated the feasibility of the proposed link unit equipped with a vectorable dual-rotor for different robot configurations. Although a similar articulated aerial robot proposed in [27] also achieved a snake-like motion in midair, the related link unit was composed of eight rotors, which led a larger size than ours. Therefore, our proposed link unit has a significant advantage in maneuvering in tight terrain that was verified in the experiment of Figure 11.

Conclusions
In this paper, we first presented a generalized design to achieve the configuration of a snake-like aerial robot where the dual-rotor vectoring apparatus is the key feature. Then, a generalized modeling method for the articulated aerial robot based on an approximated model was presented. According to the number of rotor vectoring apparatuses, underactuated and fully actuated models were further derived, which resulted in two different allocation strategies. Furthermore, a generalized control framework was developed to support both under-actuated and fully actuated models. Finally, two different platforms were built to perform experiments that involved the joint motion in midair and demonstrated the feasibility of our proposed methods.
An important improvement in the control method is the switch between the dual-rotor and virtual-single-rotor modes, as presented in Section 3. Given the switch can cause a discontinuous change in the rotor's thrust, the flight may become unstable after the transition. Then, it is necessary to develop an interpolation approach to switch between these two modes smoothly. Another key issue that remains in this work is the feasibility of our design, modeling, and control methods for configuration with more links (e.g., >4). A possible problem with a large link number is the elastic vibration due to the lack of bending and torsional rigidity of the link structure. We will utilize the redundancy of the control input to suppress the elastic vibration according to the additional control method presented in [35]. Furthermore, the self-localization is another crucial challenge for this articulated structure, and sensor fusion using the multimodal sensors distributed in each link can be developed in the future. Last but not least, more evaluations of our articulated aerial platform in the field will be performed in future work to demonstrate the advantages of snake-like structure and motion in the aerial domain. Aside from the snake-like squeezing motion for inspection in the tight and elevated terrains, the whole body grasping that imitates snake coiling can be investigated to achieve large object transportation.