Design and Motion Analysis of a Soft Modular Robot for Diverse Environments

Abstract

This article introduces the design and development of a modular soft robot capable of performing multiple movement modes. The core unit module features a four-chamber soft structure, separated by a cross-shaped thin plate. By selectively applying pneumatic pressure to different chambers and changing connector configurations, the robot achieves diverse modular configurations and movement modes, enabling it to adapt to various environments. To address the challenges posed by the material’s nonlinear behavior and its infinite degrees of freedom, a three-dimensional spatial mathematical modeling approach is proposed. This method, grounded in classical plate theory and the chained composite model, establishes a static model for the soft robot’s spatial bending motion with constant curvature. In addition, a single-controller framework based on a central pattern generator is developed to facilitate the generation of multiple movement gaits. By tuning parameters such as oscillator phase, frequency, load factor, and amplitude, the controller can generate a wide range of movement patterns. To validate the proposed theoretical and experimental models, we developed a pneumatic control platform that demonstrated the robot’s multimodal locomotion capabilities through systematic testing in terrains with varying complexity.

Introduction

The development of soft robots is intricately linked to advancements in manufacturing materials and structural design. With the rapid progress in materials science, particularly the emergence of flexible and adaptive materials, soft robots are now capable of operating in a broader range of environments and applications.^1–4 The structural design of soft robots plays a pivotal role in determining their mobility and functionality, as evidenced by various bioinspired prototypes. Examples include inchworm-like soft robots,^5,6 worm-inspired crawling robots,⁷ leech-mimicking soft robots,⁸ octopus-inspired soft robots,⁹ flagellum-inspired soft gripper,¹⁰ snake-scaled robots,¹¹ frog-inspired swimming robots,¹² starfish-inspired robots, and wall-climbing robots.¹³ These designs often incorporate flexible biomimetic structures that replicate the complex movements observed in biological organisms, enabling them to perform intricate and adaptive motions.

Flexible materials and bioinspired design principles form the foundation for the development of soft robotic systems. While these elements enable the realization of soft robotic bodies, the diversity of functionalities and breakthroughs in performance largely depend on innovative configuration designs. Traditional soft robot structures, although flexible, are often limited in functionality, prompting researchers to explore the concept of soft modular robotics.¹⁴ Soft modular robots consist of modular units and connectors, enabling changes in the structure of the modular units and the design of connectors, thereby facilitating diverse configurations and enhanced functionality.¹⁵ For instance, Xu et al.¹⁶ proposed a high-performance, inchworm-inspired modular robotic crawler utilizing fluid prestressed composite (FPC) actuators, incorporating two types of one-way wheel assemblies to adapt to diverse environments. Hu et al.¹⁷ developed a flexible modular climbing robot (Smcbot), consisting of a soft body module and three types of electromagnetic foot connectors, enabling the formation of multiple robot configurations. Zhang et al.¹⁸ introduced a soft-climbing robot combining modular pneumatic artificial muscles with two negative pressure suckers. Zhao et al.¹⁹ designed a motor rope-driven soft modular system with two types of connectors, achieving configurations such as a handgrasp and a starfish shape. MGB et al.²⁰ proposed soft building blocks inspired by elephant trunks, featuring three types of connection structures to realize configurations such as a soft hand, elephant trunk, and quadruped. Zou et al.^21,22 designed programmable soft modular units, each with unique logic functions, which can be combined to form diverse logic circuits for controlling soft hand gestures and soft turtle movements. Liu et al.²³ designed a pentagonal soft modular unit and developed a connector with a slot structure to facilitate module connections. Jin et al.²⁴ introduced mortise-and-tenon connectors, enabling modular units to form a soft gripper. These advancements in the structural design of soft modular units and connectors have significantly enriched the morphological and functional diversity of soft robots. However, the pursuit of multiple configurations has led to increasingly complex connector designs. Simplifying connector structures while ensuring seamless connections between modular units and enabling a wide range of biomimetic motions remains a critical challenge in the field.

Soft modular robots are constructed by connecting multiple soft modular units, each possessing independent movement capabilities, through connectors. The inherent nonlinearities of these soft units present significant challenges in developing accurate kinematic models for soft robots. To address these challenges, researchers have explored various approaches, including simplified modeling, the Finite Arc Method (FAM), finite element methods, bilayer plate theory, and the Chain Composite Model (CCM).^25,26 For instance, Zhang et al.²⁷ developed a simplified bending beam model to describe the continuous motion of worm-like soft robots by controlling their initial curvature. Chillara et al.²⁸ proposed an analytical model for multifunctional laminated composite materials, which consist of a mechanically prestressed layer, a fluid layer, and a constraint layer. Sayadi et al.²⁹ introduced and validated a novel mechanical model for tendon-actuated soft-bodied robots based on FAM, enabling the derivation of control equations by parameterizing deformations. Zhang et al.³⁰ established a kinematic model using bilayer plate theory, correlating the overall deformation of flexible robots with external environmental factors. In addition, Zhang et al.³¹ proposed a mathematical modeling method for chained soft robots with capsule structures, establishing the relationship between gas pressure in the cavity and the midplane deformation of soft-bodied robots. Zhou et al.³² further advanced the field by combining the CCM to model large deflections in FPC actuators. These methods have significantly contributed to the kinematic modeling of soft robots in planar motion. However, to achieve a broader range of biomimetic motion modes, it is crucial to extend these modeling approaches to three-dimensional (3D) space. While some researchers have applied the Denavit–Hartenberg (D-H) parameter method to develop 3D kinematic models for continuum robots,³³ continuum manipulators,³⁴ soft grippers,³⁵ and soft actuators,³⁶ these methods are not directly applicable to modular soft robots. The unique configuration of modular soft robots, which consist of interconnected soft units with varying degrees, necessitates specialized modeling techniques. Consequently, deriving a comprehensive kinematic model for modular soft robots in 3D space remains an open challenge, requiring further research and innovation.

Soft modular robots are composed of interconnected soft units, each equipped with two or more air chambers. Coordinating the control of these internal chambers to generate rhythmic biomimetic motions is a prominent area of research. Many researchers have explored gait control laws for these chambers through simulations and experiments. For example, Lipson et al.^37,38 utilized evolutionary algorithms in simulations to adjust expansion and contraction parameters, deriving motion control rates for soft modular units and enabling biomimetic motion. Sui et al.³⁹ enhanced the evolutionary algorithm by integrating vision, achieving obstacle avoidance motion control for modular soft robots. Peng et al.⁴⁰ experimentally implemented inchworm-inspired motion in modular soft robots and investigated the influence of control parameters on motion performance. Li et al.⁴¹ developed a PID control model for soft modular robots and, by incorporating a sensor system, achieved motion control for configurations such as quadruped, gripper, and chair. In addition, Lee et al.⁴² analyzed the deformation characteristics of soft modular units through experiments and simulations, discussing motion performance based on experimental results. While these methods successfully achieve motion control for soft robots or individual soft modular units, they are not well-suited for multimode motion control of modular soft robots. In contrast, central pattern generator (CPG) method offers the capability to control multimode movements and has been widely applied to earthworm-like robots⁴³ and snake robots⁴⁴ to generate rhythmic biomimetic motions. However, CPG output signals are typically used to control motors in modular robots, limiting their direct applicability to pneumatic soft modular robots. Leveraging the advantages of CPG control to achieve multimode motion control in pneumatic soft modular robots remains an unresolved challenge that requires further exploration.

Inspired by the continuous crawling motion of inchworms and recent advancements in pneumatic soft robotics, this article introduces a bionic-motion soft modular robot (B-SMR). The B-SMR is designed with soft modular units and a simplified connector structure, featuring a slot base and hook–claw mechanism. This design reduces structural complexity while enabling the reconfiguration of the soft modular robot. We propose a kinematic model for chained, worm-like soft modular robots to establish the relationship between the morphology of bionic motions and the driving pressure in 3D space. The accuracy of the kinematic model is validated by comparing theoretical predictions with experimental results for each motion mode. Furthermore, a multimodal motion control method based on CPG is developed, integrating the kinematic model to optimize control parameters for each bionic motion. The effectiveness of the control model is demonstrated through bionic motions in diverse environments and comparative experiments. These contributions provide significant insights into the design of simplified connector structures, the derivation of theoretical models, and the development of control strategies for soft modular robots. The B-SMR, capable of performing various bionic motions tailored to different environments, holds promise for adaptive motion pattern adjustments in complex and dynamic settings. In the future, this robot could be deployed in challenging environments such as wilderness or underwater, where it can collect information and perform tasks with enhanced versatility and adaptability.

Design and Fabrication of the B-SMR

Design of SMR

The SMR integrates two key components as follows: the soft modular unit and the connection unit. The soft modular unit contains two airbags, each divided into two chambers, which enable independent control of four air chambers through high-pressure actuation. This configuration facilitates various bionic motions. The connection unit comprises two parts as follows: a slot base (Part A) and a hook–claw (Part B). The slot base features four claw grooves distributed at 90-degree intervals, whereas the hook–claw has two claws positioned at 180-degree intervals. This arrangement allows two distinct connection modes, supporting flexible assembly and reconfiguration (Fig. 1B). Following the design phase, the fabrication stage uses 3D printing and silicone molding technologies to manufacture the soft modules and connectors.

FIG. 1.

Design components, fabrication, and bionic motions: (A) Design a modular soft robot capable of mimicking the crawling motion of inchworms. (B) Design an illustration of a soft modular robot composed of soft modular units and connection units. The soft modular unit features a symmetrical capsule structure, whereas the connection unit uses slot and hook–claw connection design. (C) Fabrication of the soft modular unit and the connection unit, and the assembly of the design components into two types of the soft modular robots (Type A and Type B). (D) Analysis of the bionic motions of the modular soft robots. The soft robots are capable of swimming, crawling, turning, and wriggling through actuated airbag. (E) More morphological configurations, such as those of cordylus cataphractus, snakes, earthworms, and ribbon eels.

Fabrication of SMR

The fabrication of the SMR utilized AB liquid silicone (Shore hardness 20) for the soft modular unit. The molds for both soft and connection units were produced using resin-based 3D printing technology. To prepare the silicone, equal volumes of components A and B were combined in a beaker and mixed at 1500 rpm for 30 s with a mechanical mixer. The homogenized mixture was then poured into three mold types as follows: airbag molds, substrate molds, and connection unit base molds. For the connection unit assembly, Parts A and B were prepositioned in their respective base molds before silicone pouring.

To eliminate air bubbles, all filled molds were placed in a vacuum chamber to create a negative pressure environment. The degassed molds were then cured in a dryer at 50°C for 3 h. After curing, the components, including airbags, substrates, and connection parts, were assembled using flexible adhesive (Smooth-On). The soft modular robot had two configuration modes achieved through its connection parts (Type A and Type B) (Fig. 1C). By actuating different configurations and internal chambers, it was capable of generating various biomimetic motions.

Bionic-inspired motion analysis of the soft robot

The soft robot mimics the motions of nature animals by selectively actuating different air chambers, whereas the connection units allow for two types of configurations. By combining different numbers of modules and using two connection methods, the soft modular robot can form a greater variety of morphological configurations. Figure 1D demonstrates that two soft modular robots can perform four distinct bionic motions, whereas Figure 1E presents a broader range of potential morphological configurations, including the ring-shaped cordylus cataphractus configuration, snake configuration, earthworm configuration, and ribbon eel configuration. This article specifically studies the bionic motions formed by two modules, where the robot can perform four bionic-inspired motions as follows: crawling, turning, wriggling, and swimming.

The soft robot configuration Type-A (with two modular units arranged vertically) is capable of crawling, wriggling, and turning motions. Crawling motion is generated by actuating the upper airbag of Unit A and two auxiliary airbags of Unit B. Turning motion is generated by simultaneously actuating the rear airbag of Unit A and the lower airbag of Unit B. Wriggling motion is generated by alternately actuating the front and rear airbags of Unit B. The soft robot configuration Type-B (with two modular units arranged horizontally) is capable of swimming motion in water. Swimming motion is generated by simultaneously actuating the upper airbags of Unit A and Unit B (Fig. 1D). To achieve these biomimetic motions, the robot’s kinematic model must first be established, followed by the development of a multimode motion control system based on deformation analysis.

Deformation Modeling of the B-SMR

The B-SMR is capable of performing various bionic motions, including crawling, wriggling, turning, and swimming. A deformation model for the B-SMR is proposed based on the piecewise constant curvature (PCC) method and the two-layer plate theory. By integrating the bionic motions with the D-H coordinate transformation method, the B-SMR centerlines for four motions in 3D space are obtained. Consequently, the relationship between the B-SMR’s bending state and air pressure is established.

Specific parameters of the B-SMR

Figure 2A illustrates the coordinate system (O-XYZ) and the physical dimensions of the B-SMR. The O-XYZ system is defined with its origin at point O, located at the center of the end plane of the soft modular unit. The XY-plane and XZ-plane represent the center horizontal plane and center vertical plane of the B-SMR, respectively. Table 1 provides detailed parameters of the soft modular unit and connection unit.

FIG. 2.

Deformation modeling of the B-SMR: (A) Parameters of the B-SMR and the pneumatic experimental platform. (B) The theoretical value curves for segmentation at 4, 7, 8, 9, and 16, and comparison with the experimental curve at a pressure of 33 KPa. (C) The principles (based on virtual rod theory) for rotation from 2D plane to 3D space. (D) Solution method for calculating bionic motion curves in 3D space.

Table 1.

Description of Specific Parameters

Explanation	Parameters	Magnitude
Total length (the soft unit)	L	91 mm
Total width (the soft unit)	C	40 mm
The width of the fluid channel	C ₁	36 mm
Total height of chamber	H	19 mm
Total height of airbag	H ₁	16 mm
The height of channel bottom	H₂	3 mm
Thickness of airbag	t ₁	2 mm
Thickness of substrate	t ₂	2.5 mm
The width of the internal rib plate	C_r	2 mm
Total length (the connection unit)	l	16 mm
The Young’s modulus of constraining layer	E ₁ E ₂	1200 MPa
The Young’s modulus of fluidic layer	E ₁ E ₂	1.2 MPa
The Poisson’s ratio	$μ_{12}$ $μ_{21}$	0.48

Strain energy of the i-th segment

The B-SMR is constructed from silicone material and performs bending motions when driven by high-pressure air. Consequently, the mathematical model of the B-SMR must account for large-deflection nonlinear deformation. The static response of the large-deflection modeling CCM is determined by dividing the actuator into multiple elements and formulating mathematical models for each element using the Small Rotation Model.³⁰

Figure 2B illustrates the discretization of the centerline of the soft modular unit into multiple segments. The strain energy model of the soft unit is based on the Lagrange strain formulation derived from classical laminated plate theory, material mechanics, and geometric nonlinearity.

By integrating the PCC principle with the CCM method, the centerline of the soft modular unit is divided into discrete segments. The global coordinate system is defined as O-XZ, whereas each segment i has its own local coordinate system O _i -X _i Y _i Z _i , with the origin positioned at the start of each segment. The undeformed length of the element before any deflection underload is denoted as L_i. The rotation angle at the free end of the i-th segment in the local coordinate system is denoted as $β_{i}$ , whereas the rotation angle at the point O _i in the global coordinate system is denoted as $α_{i}$ (Fig. 2B). Therefore, the relationship between the angles of the centerline is expressed as follows: $α_{i} = \sum_{k = 1}^{i - 1} β_{k} (i \in [2, ..., N])$ (1)

Here, the rotation angle of the first segment $α_{1}$ is set to zero.

The static response of the soft modular unit’s midplane rotation is expressed as a polynomial function with unknown coefficients, as shown below. $u_{i} = b_{i 1} x + b_{i 2} x^{3}, v_{i} = c_{i} y, w_{i} = d_{i} x^{2} (x \in [0, L_{i}], y \in [- C / 2, C / 2], i \in [2, ..., N])$ (2)

Here, $b_{i 1}$ , $b_{i 2}$ , $c_{i}$ , and $d_{i}$ are polynomial coefficients. $u_{i}$ , $v_{i}$ , and $w_{i}$ represent the elongations of the i-th segment along x _i , y _i , and z _i axes, respectively. L_i represents the virtual center length of the i-th segment, $L_{i} = L / N$ . The elongation of the first segment is assumed as zero.

The deformation of the soft unit is predominantly characterized by axial elongation and planar bending. Consequently, the effects of transverse shear stresses and transverse normal stress are disregarded. Furthermore, it can be assumed that the soft modular unit bending plate is defined as the fluid layer (inflating layer), and the other plates are defined as the constraint layer (uninflated layer) (Fig. 2B). The strain energy of the i-th segment of the soft modular unit is formulated as follows: $ϕ_{(i)} = ϕ_{C L (i)} + ϕ_{F L (i)}$ (3)

Here, $ϕ_{C L (i)}$ and $ϕ_{F L (i)}$ represent the strain energy of the constraint layer and the fluid layer in the i-th segment, respectively (Fig. 2B). The soft modular unit can be assumed to be a thin plate.

The strain energies of the constraint layer and the fluid layer in the i-th segment are expressed as follows: $ϕ_{C L (i)} = ϕ_{CLcha (i)} + ϕ_{CLrib (i)}$ (4) $ϕ_{F L (i)} = ϕ_{FLwhole (i)} - ϕ_{FLcha (i)} + ϕ_{FLrib (i)}$ (5)

Here, $ϕ_{CLcha (i)}$ and $ϕ_{CLrib (i)}$ represent the strain energy of the constraint layer for the air channel volume and the ribs volume, respectively. $ϕ_{FLwhole (i)}$ , $ϕ_{FLcha (i)}$ , and $ϕ_{FLrib}$ represent the strain energy of the fluid layer for the total volume, the fluid channel volume, and the rib volume, respectively.

By integrating the two-layer plate method and mechanics of materials, the strain energy of each layer is obtained as follows: $ϕ_{CLcha (i)} = \int_{- t_{2}}^{t_{2}} \int_{- C / 2}^{C / 2} \int_{0}^{L_{i}} (\frac{1}{2} Q_{11}^{(C L)} ε_{i x}^{2} + Q_{12}^{(C L)} ε_{i x} ε_{i y} + \frac{1}{2} Q_{22}^{(C L)} ε_{i y}^{2}) d x d y d z$ (6) $ϕ_{CLrib (i)} = \int_{- H / 2}^{- t_{2}} \int_{- C_{r} / 2}^{C_{r} / 2} \int_{0}^{L_{i}} (\frac{1}{2} Q_{11}^{(C L)} ε_{i x}^{2} + Q_{12}^{(C L)} ε_{i x} ε_{i y} + \frac{1}{2} Q_{22}^{(C L)} ε_{i y}^{2}) d x d y d z$ (7) $ϕ_{FLwhole (i)} = \int_{t_{2}}^{H_{2} + H_{1} / 2} \int_{- C / 2}^{C / 2} \int_{0}^{L_{i}} (\frac{1}{2} Q_{11}^{(F L)} ε_{i x}^{2} + Q_{12}^{(F L)} ε_{i x} ε_{i y} + \frac{1}{2} Q_{22}^{(F L)} ε_{i y}^{2}) d x d y d z$ (8) $ϕ_{FLcha (i)} = \int_{t_{2}}^{H_{2} + H_{1} / 2} \int_{- C_{1} / 2}^{C_{1} / 2} \int_{0}^{L^{'}_{i}} (\frac{1}{2} Q_{11}^{(F L)} ε_{i x}^{2} + Q_{12}^{(F L)} ε_{i x} ε_{i y} + \frac{1}{2} Q_{22}^{(F L)} ε_{i y}^{2}) d x d y d z$ (9) $ϕ_{F L rib (i)} = \int_{t_{2}}^{H_{2} + H_{1} / 2} \int_{- C_{r} / 2}^{C_{r} / 2} \int_{0}^{L^{'}_{i}} (\frac{1}{2} Q_{11}^{(F L)} ε_{i x}^{2} + Q_{12}^{(F L)} ε_{i x} ε_{i y} + \frac{1}{2} Q_{22}^{(F L)} ε_{i y}^{2}) d x d y d z$ (10)

Here, $L_{i}^{'}$ represents the length of the fluid channel in the i-th segment, $L_{i}^{'} = (L - 2 r_{1}) / N$ . $C_{r}$ denotes the width of the internal rib plate of the soft modular unit. $Q$ is the stiffness matrix of the material. $ε_{i x}$ and $ε_{i y}$ denote the tensile strain in the x-direction and y-direction, respectively. $γ_{ixy}$ represents the shear strain on the o_i-x_iy_i plane.

The stiffness parameters of the material are determined as follows: $Q_{11} = \frac{E_{1}}{1 - μ_{12} μ_{21}}, Q_{12} = \frac{μ_{12} E_{2}}{1 - μ_{12} μ_{21}}, Q_{22} = \frac{E_{2}}{1 - μ_{12} μ_{21}}$ (11)

Here, E₁ and E₂ represent the Young’s modulus of the material, whereas $μ_{12}$ and $μ_{21}$ denote the Poisson’s ratio of the material. The Young’s modulus of the fluidic layer differs from that of the constraining layer. Material parameters are shown in Table 1.

The tensile strain and the shear strain are calculated as follows: ${\begin{cases} ε_{i x} = \frac{\partial u_{i}}{\partial x} + \frac{1}{2} {(\frac{\partial w_{i}}{\partial x})}^{2} - z (\frac{\partial^{2} w_{i}}{\partial x^{2}}) \\ ε_{i y} = \frac{\partial v_{i}}{\partial y} + \frac{1}{2} {(\frac{\partial w_{i}}{\partial y})}^{2} - z (\frac{\partial^{2} w_{i}}{\partial y^{2}}) \\ γ_{ixy} = \frac{\partial v_{i}}{\partial x} + \frac{\partial u_{i}}{\partial y} + \frac{\partial w_{i}}{\partial x} \frac{\partial w_{i}}{\partial y} - 2 z (\frac{\partial^{2} w_{i}}{\partial x \partial y}) \end{cases} (i \in [2, ...., N])$ (12)

Pressure work

The bionic motions of the B-SMR are realized by actuating the soft robot with air pressure. The inflation of the airbag induces a volumetric change in the soft modular robot, based on the principle of conservation of energy. $\frac{\partial (\sum_{1}^{N} ϕ_{(i)}) - \partial (\sum_{1}^{N} W_{P (i)}) - \partial W_{F}}{\partial k_{i}} = 0$ (13)

Here, $W_{P}$ and $W_{F}$ present the work done on the i-th segment by air pressure and external load, respectively. By applying the statics principles, $W_{F}$ is set to zero. Furthermore, the B-SMR operates at room temperature, with the thermodynamic deformation neglected. Consequently, the work done by the pressure is fully converted into strain energy during actuation.

The work done by air pressure on the i-th segment is expressed as follows: $W_{P (i)} = \frac{P_{0} V_{i} - P_{1} V_{i}^{'}}{κ - 1}$ (14)

Here, $P_{0}$ and $P_{1}$ denote the initial and final pressure values, respectively. $V_{i}$ and $V_{i}^{'}$ represent the initial and final volumes of the soft modular units, respectively. $κ$ denotes the ratio of specific heats, $κ = 1.4$ .

The final volume of the soft modular robot is expressed as follows: $V_{i}^{'} = \int_{V_{i}} (1 + ε_{i x}) (1 + ε_{i y}) d V_{i} (i \in [2, ..., N])$ (15)

Numerical solution

The theoretical values were obtained by implementing the equations in MATLAB using the variational Rayleigh–Ritz method and the Newton–Raphson method to derive the centerline of the soft modular unit. The experimental values were obtained by inflating the soft modular unit with air on the experimental platform. The pneumatics experimental platform consists of an air pump, power supply, relay, solenoid valves, microcontroller, pressure sensors, and PC, as shown in Figure 2A. Figure 2B shows both the theoretical and experimental curves with the air pressure set at 33 kPa. The theoretical curves were obtained by solving with segmentation numbers of 4, 7, 8, 9, and 16, whereas the experimental curve was determined through experimental measurements and the point plotting method. As the number of segments increases, the theoretical results approach the experimental results. However, increasing the number of segments does not always result in better alignment with the experimental data. As the number of segments increases, the curve’s convergence strengthens. When the number of segments is less than eight, the theoretical values continuously approach the experimental results. Conversely, when the number of segments exceeds eight, the theoretical values increasingly deviate from the experimental results. Therefore, for the soft modular robot’s bionic-motion theoretical curves, the number of segments is set to eight.

Theoretical curves for bionic motions from two-dimensional plane to 3D space

Using the aforementioned modeling method, the static model of the soft modular unit was successfully established on a two-dimensional (2D) plane. Given that each modular unit shares identical dimensions, material properties, and actuation mechanisms, the derived CCM static model on the 2D plane is generalizable. However, the modular soft robot incorporates two connection configurations (Type A and Type B), requiring an extension of its bionic motion analysis into 3D space. As shown in Figure 2C, the centerlines of the bionic motions were extended from 2D to 3D by integrating the D-H parameters method and the virtual rod principle. The yellow line represents the centerline in the 2D O-XZ plane, whereas the final bionic motion is depicted within the $O - X' Y' Z'$ 3D coordinate system. Using the virtual rod principle, the virtual rod connects the origin of the coordinate system to points on the theoretical curve. By applying the virtual rod principle, the virtual rod links the origin of the coordinate system to points on the theoretical curve. Utilizing the D-H parameters method, the virtual rod first rotates around the $Y^{'}$ -axis, followed by a rotation around either the $Z^{'}$ -axis or $X^{'}$ -axis. Finally, the theoretical results are visualized in the 3D coordinate system.

The 3D coordinate representation of the soft modular robot system is defined in accordance with the D-H convention as follows: $(\begin{matrix} X_{i}^{'} \\ Y_{i}^{'} \\ Z_{i}^{'} \\ 1 \end{matrix}) = (\begin{matrix} R_{x} (δ) & 0 \\ 0 & 1 \end{matrix}) (\begin{matrix} R_{z} (φ) & 0 \\ 0 & 1 \end{matrix}) (\begin{matrix} R_{y} (θ_{i}) & 0 \\ 0 & 1 \end{matrix}) (\begin{matrix} 0 \\ 0 \\ ξ_{i} \\ 1 \end{matrix})$ (16)

The rotation matrices are defined as follows: $R_{y} (θ_{i}) = (\begin{matrix} \cos θ_{i} & 0 & \sin θ_{i} \\ 0 & 1 & 0 \\ - \sin θ_{i} & 0 & \cos θ_{i} \end{matrix})$ (17) $R_{z} (φ) = (\begin{matrix} \cos φ & - \sin φ & 0 \\ \sin φ & \cos φ & 0 \\ 0 & 0 & 1 \end{matrix})$ (18) $R_{x} (δ) = (\begin{matrix} 1 & 0 & 0 \\ 0 & \cos δ & - \sin δ \\ 0 & \sin δ & \cos δ \end{matrix})$ (19)

Here, $X_{i}^{'}$ , $Y_{i}^{'}$ , and $Z_{i}^{'}$ represent the coordinates of point i in the $O - X' Y' Z'$ coordinate system. $ξ_{i}$ denotes the length of the virtual rod of point i in the O-XZ plane. $θ_{i}$ represents the rotation angle of the point i around the $Y^{'}$ -axis, whereas $φ$ and $δ$ denote the rotation angles around the $Z^{'}$ -axis and $X'$ -axis, respectively, within the $O - X' Y' Z'$ coordinate system.

The rotation angle $θ_{i}$ is obtained as follows: $θ_{i} = \tan (\frac{Z_{i}}{X_{i}}) (i \in [2, ..., N])$ (20)

Here, $X_{i}$ and $Z_{i}$ represent the coordinates of point i in the O-XZ plane.

Crawling motion. Crawling motion is suitable for moving forward on flat surfaces. The B-SMR performs crawling motion when the upper airbag of unit A is inflated and the airbags of unit B are inflated with assistance. To facilitate the description and analysis of this motion, the origin of the coordinate system is established at the start of unit A, as illustrated in Figure 2D(i). In this configuration, the initial line OA, which originally lies within the O-XZ plane, can be transformed into the $O - Y' Z'$ plane. The coordinates are defined as follows: $(\begin{matrix} X_{i}^{'} \\ Y_{i}^{'} \\ Z_{i}^{'} \end{matrix}) = (\begin{matrix} 0 \\ ξ_{i} \sin θ_{i} \\ ξ_{i} \cos θ \end{matrix}), θ_{i} \leq 0$ (21)

Here, $X_{i}^{'}$ , $Y_{i}^{'}$ , and $Z_{i}^{'}$ represent coordinates of the point i on line OA. $φ = π / 2$ and $δ = 0$ .

The lengths of the remaining centerline segments can be derived equivalently from the dimensions of the connector and Unit B. Consequently, the coordinates of points B and C are defined as follows: $B_{X' Y' Z'} (X_{A}^{'}, Y_{A}^{'} + l \cos (α_{N} - \frac{π}{2}), Z_{A}^{'} + l \sin (α_{N} - \frac{π}{2}))$ (22) $C_{X' Y' Z'} (X_{B}^{'}, Y_{B}^{'} + L \cos (α_{N} - \frac{π}{2}), Z_{B}^{'} + L \sin (α_{N} - \frac{π}{2}))$ (23)

Wriggling motion. Wriggling motion is suitable for movement on soil-covered ground, resembling the serpentine locomotion of snakes. This motion is generated by alternately actuating the front and rear airbags of Unit B while maintaining continuous actuation of the upper airbag of Unit A. The origin of the coordinate system is established at the starting point of Unit A, as shown in Figure 1D(ii). For instance, by actuating the rear airbag of Unit B and the upper airbag of Unit A, the robot achieves this wriggling motion. The calculation methods for lines OA and point B are consistent with those used in crawling motion and can be derived from Eq. (21) and Eq. (22), respectively. The coordinates of line BC are defined as follows: $(\begin{matrix} X_{j}^{'} \\ Y_{j}^{'} \\ Z_{j}^{'} \end{matrix}) = (\begin{matrix} ξ_{i} \sin θ_{i} \\ - ξ_{i} \cos θ_{i} \sin δ \\ ξ_{i} \cos θ_{i} \cos δ \end{matrix}) + (\begin{matrix} X_{B}^{'} \\ Y_{B}^{'} \\ Z_{B}^{'} \end{matrix})$ (24)

Here, $X_{j}^{'}$ , $Y_{j}^{'}$ , and $Z_{j}^{'}$ present coordinates of the point j on line ${BC}^{'}$ . $φ = 0$ and $δ = - α_{N}$ .

Turning motion. The turning motion is generated by actuating either the upper or lower airbags of Unit B while maintaining continuous actuation of the rear airbag of Unit A. The origin of the coordinate system is established at the starting point of Unit A, as shown in Figure 1D(iii). For instance, actuating the lower airbag of Unit B and the rear airbag of Unit A results in an anticlockwise rotation around the $Z^{'}$ -axis. The first line OA can be transformed from the O-XZ plane to the $O - Y' Z'$ plane. The coordinates are expressed as follows: $(\begin{matrix} X_{i}^{'} \\ Y_{i}^{'} \\ Z_{i}^{'} \end{matrix}) = (\begin{matrix} ξ_{i} \cos θ_{i} \\ ξ_{i} \sin θ_{i} \\ 0 \end{matrix})$ (25)

Here, $X_{i}^{'}$ , $Y_{i}^{'}$ , and $Z_{i}^{'}$ present coordinates of the point i on line OA. $φ = 0$ and $δ = - π / 2$ .

The line AB can be equivalently determined by the length of the connector. The coordinates of point B are defined as follows: $B_{X' Y' Z'} (X_{A}^{'} + l \sin (α_{N} - \frac{π}{2}), Y_{A}^{'} + l \cos (α_{N} - \frac{π}{2}), Z_{A}^{'})$ (26)

The line BC can be transformed from the O-XZ plane to the $O - Y' Z'$ plane by first rotating around the $Y^{'}$ -axis, then around the $X^{'}$ -axis, and finally around the $Z^{'}$ -axis. The coordinates of line BC are defined as follows: $(\begin{matrix} X_{j}^{'} \\ Y_{j}^{'} \\ Z_{j}^{'} \end{matrix}) = (\begin{matrix} - ξ_{i} \cos θ_{i} \cos φ \\ ξ_{i} \cos θ_{i} \sin φ \\ ξ_{i} \sin θ_{i} \end{matrix}) + (\begin{matrix} X_{B}^{'} \\ Y_{B}^{'} \\ Z_{B}^{'} \end{matrix})$ (27)

Here, $X_{j}^{'}$ , $Y_{j}^{'}$ , and $Z_{j}^{'}$ represent coordinates of the point j on line BC. $θ = - π / 2$ and $φ = - α_{N}$ .

Swimming motion. The swimming motion is ideally suited for underwater movement. To maintain balance, the driving pressure values of the two soft units are set to be equal. The origin of the coordinate system is established at the midpoint of the connector, as shown in Figure 2D(iv). The centerline of the B-SMR exhibits symmetry with respect to the O- $X^{'} Z^{'}$ plane. Consequently, the translation process initiates with movement in the O-XZ plane, followed by sequential rotations around the $Y^{'}$ axis and the $Z^{'}$ axis. The final results are then expressed in O- $X^{'} Y^{'}$ plane. The coordinates of points B and C are $(0, - l / 2, 0)$ and $(0, l / 2, 0)$ , respectively. The coordinates of lines AB and CD are defined as follows: $(\begin{matrix} X_{j}^{'} \\ Y_{j}^{'} \\ Z_{j}^{'} \end{matrix}) = = (\begin{matrix} ξ_{i} \cos θ_{i} \cos φ \\ ξ_{i} \cos θ_{i} \sin φ \\ - ξ_{i} \sin θ_{i} \end{matrix}) + (\begin{matrix} 0 \\ k \cdot \frac{l}{2} \\ 0 \end{matrix}), {\begin{matrix} φ = \frac{π}{2}, k = - 1 line AB \\ φ = \frac{- π}{2}, k = 1 line CD \end{matrix}$ (28)

Here, $X_{j}^{'}$ , $Y_{j}^{'}$ , and $Z_{j}^{'}$ represent the point coordinates of the point j on lines AB and CD and $θ = π / 2$ .

Control Method for the B-SMR

The centerlines of bionic motions were established through kinematic modeling, a methodology that fundamentally defines the pressure-motion state correlation. Moreover, the pneumatic actuation sequence exhibits critical influence on locomotor performance. Therefore, the controller requires precise modulation of the robot’s motion modes, bending amplitude, and cycle durations. To address this, the CPG method⁴³ and Hopf oscillator⁴⁴ are introduced into the B-SMR controller. Figure 3A illustrates that each oscillator is composed of two mutually coupled neurons. The B-SMR consists of two soft modular units, each equipped with two airbags. Furthermore, each airbag requires two air tubes to control inflation and deflation. Therefore, the CPG controller consists of eight oscillators (Fig. 3B), and its governing equations are defined by the following second-order differential equation: ${\begin{cases} {\overset{\cdot}{u}}_{i} = - ω_{i} v_{i} + σ (μ - u_{i}^{2} - v_{i}^{2}) u_{i} \\ {\overset{\cdot}{v}}_{i} = ω_{i} u_{i} + σ (μ - u_{i}^{2} - v_{i}^{2}) v_{i} \end{cases} (i \in [1, ..., 8])$ (29)

FIG. 3.

CPG controller for the B-SMR: (A) The output signal of CPG network. (B) The CPG control network consists of eight oscillators (①–④ four oscillators for inflation, ⑤–⑧ four oscillators for deflation). The output signals are processed in three steps to generate the airbag control signals. (C) Airbag actuation sequence for crawling, wriggling, turning, and swimming motions.

Here, u_i and v_i are the state variables, whereas $ω_{i}$ and $μ$ represent the oscillator’s angular frequency and square of the amplitude, respectively, with the amplitude $A = \sqrt{μ}$ and $μ > 0$ .

The output signals of the CPG network are shown in Figure 3A. $σ$ is the constant that controls the convergence rate of u_i and v_i toward the limit cycle. When the state variables u_i and v_i take nonzero values, a stable limit cycle exists, defined by $\sqrt{u_{i}^{2} + v_{i}^{2}} = \sqrt{μ}$ , enabling the oscillator to produce a stable periodic oscillatory signal. If the initial values of the state variables are zero, the output remains stable at the origin. Although the output of the CPG controller stabilizes over time, it cannot be directly used to control the movement of the B-SMR. Therefore, the CPG controller is processed in three steps to enable the sequential activation of the airbags in the B-SMR. The first step introduces a phase difference, as follows: $(\begin{matrix} {\overset{\cdot}{u}}_{i} \\ {\overset{\cdot}{v}}_{i} \end{matrix}) = (\begin{matrix} σ (μ - u_{i}^{2} - v_{i}^{2}) & - ω_{i} \\ ω_{i} & σ (μ - u_{i}^{2} - v_{i}^{2}) \end{matrix}) (\begin{matrix} u_{i} \\ v_{i} \end{matrix}) + λ \sum_{j} (\begin{matrix} \cos θ_{j}^{i} & - \sin θ_{j}^{i} \\ \sin θ_{j}^{i} & \cos θ_{j}^{i} \end{matrix}) (\begin{matrix} u_{i} \\ v_{i} \end{matrix}) (i \in [1, ..., 8])$ (30)

Here, $λ$ denotes the coupling strength between oscillators, whereas $θ_{j}^{i}$ denotes the phase difference between the i^–th and j^–th oscillators. As shown in Figure 3A, the output signals of the CPG controller governing the inflation/deflation of the eight airbags exhibit phase differences, resulting in nonoverlapping trajectories. This generates a time-symmetric phase-coordinated gait. Furthermore, through parameter tuning of $λ$ , $σ$ , and $θ_{i}^{j}$ , symmetrical drive signals with uniform spatial distribution and identical amplitudes can be obtained, as demonstrated in Figure 3B.

The second step is refined to incorporate the rise and fall times, as follows: ${\begin{cases} ω_{i} = \frac{ω_{s t}}{e^{- a v_{i}} + 1} + \frac{ω_{s w}}{e^{a v_{i}} + 1} \\ ω_{s t} = \frac{1 - τ}{τ} ω_{s w} \end{cases}$ (31)

Here, $ω_{s w}$ and $ω_{s t}$ regulate the angular frequency of the rise time and fall times, respectively. $a$ represents the switching speed between $ω_{s t}$ and $ω_{s w}$ . $τ$ denotes the degree of asymmetry in the output signal, which is determined by the load factor. The parameters $ω_{s w}$ and $τ$ together influence the output signal period T, with the relationship given by $T = π τ / (1 - τ) ω_{s w} + π / ω_{s w}$ . These new gaits, characterized by asymmetrical rise and fall phases, are referred to as time-asymmetrical phase-coordinated gaits (Fig. 3B).

The third step involves discretizing the signals and converting them into control signals for the inflation/deflation of each airbag. A binary vector $w_{i} \in {0, 1}$ , used for pneumatic valve actuation in this step, is defined as follows: $w_{i} = {\begin{cases} 1, d u_{i} / d t_{i} > 0 \\ 0, else \end{cases} (i \in [1, ..., n])$ (32)

The continuous-time transformation signal is mapped into a time-domain step signal using the Heaviside step function to construct an activation function. Mathematically, the time derivative of u_i, denoted as $d u_{i} / d t_{i}$ , describes the signal’s variation, generating a binary vector in the form of a rectangular wave for gait discretization. Negative parameters are assigned a value of “0,” and positive parameters are assigned a value of “1,” as shown in Figure 3B.

In addition, a constraint is applied such that the sum of the rise time ratio $r_{rise}$ and the fall time ratio $r_{fall}$ equals 1. This ensures that the signal period remains constant while preserving the relationship between the rise time $t_{rise}$ and the fall time $t_{fall}$ , as expressed by: ${\begin{cases} r_{rise} + r_{fall} = 1 \\ t_{rise} + t_{fall} = (r_{rise} + r_{fall}) T = T \end{cases}$ (33)

Figure 3C illustrates the pneumatic actuation sequence governing biomimetic locomotion. Crucially, auxiliary chambers (depicted in light green) are integrated into the crawling and wriggling motion modules to suppress parasitic rotational dynamics. Without these chambers, the system’s kinematic behavior degrades into uncontrolled rotation or remaining stationary. Furthermore, distinct phase synchronization patterns characterize the actuation sequences for each biomimetic mode, as quantified by the CPG controller parameters listed in Table 2. (Simulation results of the CPG controller are shown in Supplementary Data).

Table 2.

Controller Parameters

Explanation	$T (s)$	$A (°)$	$v$ (velocity)
Crawling	9.2	119	$2.44 mm / s$
Wriggling (sand)	10.2	129	$0.64 mm / s$
Wriggling (soil)	10.2	128	$\begin{array}{l} 1.88 mm / s \end{array}$
Turning (clockwise)	6.8	129	$2.2 ° / s$
Turning (anticlockwise)	6.8	131	$1.84 ° / s$
Swimming	3.4	125	$17.89 mm / s$

Pressure Difference: 0.26 MPa.

Results and Discussion

Theoretical model verification

By integrating the material properties and structural parameters of the soft robot, as detailed in Table 1, an optimal segment count of N = 8 was selected for modeling the soft robot (Fig. 2B). Figure 4 illustrates the centerlines of the B-SMR biomimetic motions in 3D space, derived from the kinematic model within a 3D coordinate system. The 3D models of the B-SMR’s biomimetic motions are shown at an air pressure of 33 kPa. This pressure value was chosen because experimental results demonstrate that motion performance across all modes is optimal at this specific pressure, as further discussed in Supplementary Data.

FIG. 4.

Comparison of experimental and theoretical curves for various biomimetic motions of the B-SMR.

The theoretical curves align well with the experimental curves for the four biomimetic motions, demonstrating similar trends and confirming the effectiveness of the modeling method. Specifically, the theoretical and experimental values for the swimming motion are nearly identical, with almost complete overlap. In water, buoyancy fully counteracts gravity, and the symmetric configuration of the robot (Type B) used for swimming ensures that the theoretical model closely matches the experimental results, with an overall error of less than 1 degree. However, the curves for the other three motions exhibit some discrepancies. The maximum errors for these motions are $8^{\circ}$ , ${6.5}^{\circ}$ , and ${12.5}^{\circ}$ . These discrepancies primarily arise from the use of asymmetrical configurations (Type A) for the robot during these motions. In addition, the theoretical model, which applies the constant curvature method to process the connectors, did not account for the full impact of the connectors, which could not be fully predicted.

Motion experiments

To determine the optimal speed of the B-SMR, a driving air pressure of 33 kPa was applied across all experimental motions. Figure 5 presents the continuous waveform, the binarized signal from the CPG, and the results of the biomimetic motion experiments. The continuous waveform illustrates the movement pattern of the soft robot, whereas the binarized signal is used to control the solenoid valve, enabling rhythmic biomimetic motion.

FIG. 5.

The output signals of the CPG controller and the bionic motion experiments of the B-SMR. CPG.

Crawling motion is ideal for movement on flat surfaces. The robot achieves movement by increasing its contact area with the ground to enhance static friction, thereby propelling the soft robot forward (Supplementary Movie S1). The time duration for each crawling cycle is 9.2 s, and the average crawling speed is 2.44 mm/s, measured after crawling a distance of 27 centimeters over 12 cycles.

The turning motion includes both clockwise and counterclockwise rotations (Supplementary Movie S2). Each turning cycle, regardless of direction, lasts 6.8 s. The B-SMR completes a 180-degree turn after 12 cycles, achieving an average angular velocity of 2.2 degrees per second for both clockwise and counterclockwise rotations.

The wriggling motion is particularly effective for navigating soil- or sand-covered terrains. In the experiments, a 10 mm layer of soil and sand was used to simulate these environments (Supplementary Movie S3). Each wriggling cycle lasts 10.2 s, with the robot achieving a faster wriggling speed on soil-covered ground (1.88 mm/s) compared with sand-covered ground (0.64 mm/s).

In the swimming motion, the robot moves through water (Supplementary Movie S4). Each swimming cycle lasts 3.4 s and is repeated over 6 cycles, achieving an average speed of 17.89 mm/s. However, due to the submerged air supply pipeline, the inflation process becomes uneven. In addition, the robot encounters complex hydrodynamic forces, which cause its movement trajectory to deviate from a straight line.

Comprehensive experiment

Figure 6 demonstrates the B-SMR operating in a comprehensive experimental environment, which includes both flat ground and soil-covered terrain (Supplementary Movie S5). The experiment is divided into two stages, involving crawling, turning, and wriggling motions. In the first stage, the B-SMR moves on flat ground. It begins with three cycles of crawling motion to approach the soil-covered area, followed by three cycles of clockwise turning motion along the edge of the soil. Next, it performs three cycles of wriggling motion and concludes with another three cycles of clockwise turning motion. In the second stage, the B-SMR transitions to the soil-covered ground. It initially executes 12 cycles of wriggling motion, then repeats three cycles of counterclockwise turning motion. After completing 11 additional cycles of wriggling motion, it exits the soil-covered area. Following this, it performs three more cycles of counterclockwise turning motion and finishes with two cycles of crawling motion along the soil’s edge.

FIG. 6.

Comprehensive experiment. (A) Schematic diagram of the motion process. (B) The experiment consists of two stages as follows: moving on flat ground (a–d) and moving on soil-covered ground (e–i).

The average movement speed and rotation speed of each biomimetic motion in the comprehensive movement are generally consistent. This indicates that the robot’s mobility is stable and further validates the effectiveness of its theoretical derivation and motion control.

Conclusion

This article presents a modular soft robot capable of multimotion modes. Using classical plate theory and constant curvature chain (CCM) modeling techniques, the centerline of the soft actuator unit module in a 2D plane was derived. The proposed kinematic model maps virtual rigid rods to a spatial coordinate system and incorporates an intermediate rigid link into the mathematical framework. To streamline coordinate transformation, a local coordinate system was introduced, enabling the development of a static centerline for the constant curvature chain model in 3D space. For bionic motion control, a CPG-based neural network model incorporating a Hopf oscillator was proposed to regulate the soft robot’s movements. The model uses parameters such as load factor, frequency, and amplitude to control the robot’s various motion modes, facilitating multimode locomotion across diverse environments. Finally, the accuracy of the theoretical and control models was validated through comparisons with experimental centerline results and motion control tests.

The correctness of the theoretical model is confirmed through experimental comparison. The centerline derived from the theoretical model closely matches the experimental centerline, with maximum errors in crawling, turning, and wriggling motions being $8^{\circ}$ , ${6.5}^{\circ}$ , and ${12.5}^{\circ}$ , respectively. In swimming, the theoretical centerline aligns closely with the experimental centerline, with an error of less than 1 degree. The effectiveness of the controller is demonstrated through the robot’s bionic motion and comprehensive experiments. The speeds for crawling, turning, anticlockwise wriggling on soil-covered ground, wriggling on sand-covered ground, and swimming are 2.44 mm/s, ${2.2}^{\circ} / s$ , 1.88 mm/s, 0.64 mm/s, and 17.89 mm/s, respectively. Experimental results indicate that when the soft robot operates in challenging environments, such as sandy terrain or water, its movement speed decreases or the controllability of its motion direction is limited.

Future research could focus on improving the robot’s connection parts and control methods to support more morphological configurations. By incorporating electrification design, we can address the limitations of manually changing the robot’s configuration. The B-SMR enables a broader range of potential morphological configurations, including the ring-shaped cordylus cataphractus configuration, snake configuration, earthworm configuration, and ribbon eel configuration. Moreover, as more morphological configurations consist of many soft units, coordinating the control of the internal chambers across different modules to generate effective biomimetic rhythmic motions becomes even more challenging. Therefore, optimizing the robot’s structure and conducting research on multimode modular robot control methods in future work will enhance the robot’s capabilities for environmental exploration.

Footnotes

Authors' Contributions

Y.Z. conceived the idea and supervised the research, and revised the article. Y.L. designed prototypes, processed the data and drafted this article. D.S., L.L., and T.Z. drew and revised the figures. Z.Z., S.Z., and F.Z. conducted the numerical simulations. D.L. and Y.Z. collected and processed the data.

Author Disclosure Statement

No competing financial interests exist.

Funding Information

The study was funded and supported by National Natural Science Foundation of China (Grant no. 52025054, 52435001), Heilongjiang Province Key Research and Development Program of China (Grant no. JD2023SJ18), and State Key Laboratory of Robotics and Systems (HIT) (SKLRS-2023-KF-18, SKLRS-2025-KF-13).

Supplementary Material

References

Baines

, Patiballa

, Booth

, et al. Multi-environment robotic transitions through adaptive morphogenesis. Nature, 2022; 610(7931):283–289.

, Wong

, Shih

, et al. Bioinspired soft robots for deep-sea exploration. Nat Commun, 2023; 14(1):7097.

Zhang

, Lei

, Chen

, et al. Integrated 3D printing of flexible electroluminescent devices and soft robots. Nat Commun, 2022; 13(1):4775.

Yao

, Kim

, Yin

, et al. Multimodal soft robotic actuation and locomotion. Adv Mater, 2024; 36:2308829.

Chen

, Hu

, Zou

, et al. Design and fabrication of a multi-motion mode soft crawling robot. J Bionic Eng, 2020; 17(5):932–943.

Qin

, Liang

, Huang

, et al. A versatile soft crawling robot with rapid locomotion. Soft Robot, 2019; 6(4):455–467.

, Wagner

, Liu

, et al. Locomotion of an untethered, worm-inspired soft robot driven by a shape-memory alloy skeleton. Sci Rep, 2022; 12(1):12392.

, Liu

, Hu

, et al. Soft magnetic microrobot doped with porous silica for stability-enhanced multimodal locomotion in a nonideal environment. ACS Appl Mater Interfaces, 2022; 14(8):10856–10874.

Wehner

, Truby

, Fitzgerald

, et al. An integrated design and fabrication strategy for entirely soft, autonomous robots. Nature, 2016; 536(7617):451–455.

10.

Jiao

, Zhu

, Hu

, et al. Flagellum-inspired soft rotary motor. Adv Mater Technol, 2023; 8(18):2300476.

11.

Mei

, Yu

, Tang

, et al. A soft crawling robot with multi-modal locomotion inspired by the movement mechanism of snake scales. IEEE Robot Autom Lett, 2024; 9(9):7589–7596.

12.

Fan

, Wang

, Yu

, et al. Swimming performance of the frog-inspired soft robot. Soft Robot, 2020; 7(5):615–626.

13.

Huang

, Zhang

, Feng

, et al. A bionic starfish adsorption crawling soft robot. J Bionic Eng, 2024; 21(1):149–165.

14.

Jiao

, Zhang

, Hong

, et al. Piezo-wormbots for continuous crawling. Soft Robot, 2024; 11(2):260–269.

15.

Zhang

, Zhu

, Lin

, et al. Modular soft robotics: Modular units, connection mechanisms, and applications. Adv Intell Syst, 2020; 2:1900166.

16.

, Hu

, Xiao

, et al. Modular soft robotic crawlers based on fluidic prestressed composite actuators. J Bionic Eng, 2024; 21(2):694–706.

17.

, Dong

, Sun

. Soft modular climbing robots. IEEE Trans Robot, 2023; 39(1):399–416.

18.

Zhang

, Yang

, Yan

, et al. Inchworm inspired multimodal soft robots with crawling, climbing, and transitioning locomotion. IEEE Trans Robot, 2022; 38(3):1806–1819.

19.

Zhao

, Jiang

, She

, et al. SoftSnap: Rapid prototyping of untethered soft robots using snap-together modules. Press arXiv; 2410:19169.

20.

Atia

MGB

, Mohammad

, Gameros

, et al. Reconfigurable soft robots by building blocks. Adv Sci, 2022; 9:2203217.

21.

Jiao

, Hu

, Dong

, et al. Reprogrammable metamaterial processors for soft machines. Adv Sci, 2024; 11(11):2305501.

22.

Jiao

, Hu

, Shi

, et al. Reprogrammable, intelligent soft origami LEGO coupling actuation, computation, and sensing. Innovation (Camb), 2024; 5(1):100549.

23.

Liu

, Li

, Huang

, et al. Earthworm-inspired multimodal pneumatic continuous soft robot enhanced by winding transmission. Cyborg Bionic Syst, 2025; 6:204; doi: 10.34133/cbsystems.0204

24.

Jin

, Wang

, Xiong

, et al. Modular soft robot with origami skin for versatile applications. Soft Robot, 2023; 10(4):785–796.

25.

Mei

, Zhao

, Tang

, et al. A review of soft robot modeling and control. Robot, 2024; 46:234–256.

26.

Yasa

, Toshimitsu

, Michelis

, et al. An overview of soft robotics. Annu Rev Control Robot Auton Syst, 2023; 6(1):1–29.

27.

Zhang

, Xu

, Liu

. Dynamic modeling for continuous locomotion of bionic soft worm-like robot. Int J Non Linear Mech, 2024; 162:104702.

28.

Chillara

VSC

, Headings

, Dapino

. Multifunctional composites with intrinsic pressure actuation and prestress for morphing structures. Compos Struct, 2016; 157:265–274.

29.

Sayadi

, Cecere

, Hooshiar

. Finite arc method: Fast-solution extended piecewise constant curvature model of soft robots with large variable curvature deformations. Robot Rep, 2024; 2(1):49–64.

30.

Zhang

, Wang

, Zhou

, et al. A kinematic model for soft robot based on double laminate plate theory. In: IEEE International Conference on Robotics and Biomimetics. IEEE, Kuala Lumpur: Malaysia; 2018; pp.600–605.

31.

Zhang

, Huang

, Li

, et al. Multi-modal bionic motion analysis of a cpg-controlled pneumatic soft robot. J Bionic Eng, 2024; 21(5):2247–2257.

32.

Zhou

, Headings

, Dapino

. Modeling of fluidic prestressed composite actuators with application to soft robotic grippers. IEEE Trans Robot, 2022; 38(4):2166–2178.

33.

Zhang

, Fang

, Xiang

, et al. A survey on design, actuation, modeling, and control of continuum robot. Cyborg Bionic Syst, 2022; 2022:9754697.

34.

Wang

, Frazelle

, Wagner

, et al. Dynamic control of multisection three-dimensional continuum manipulators based on virtual discrete-jointed robot models. IEEE/ASME Trans Mechatron, 2021; 26(2):777–788.

35.

Lin

, Shahabi

, Yang

, et al. Parametrically modeled DH table for soft robot kinematics: case study for a soft gripper. In: the 15th IFToMM World Congress on Mechanism and Machine Science. Springer: Krakow, Poland; 2019;pp.617–625.

36.

Zhang

, Wei

, Bi

. Design and experimental research of the hybrid-driven soft robot. IR, 2023; 50(4):648–658.

37.

Cheney

, Lipson

. Topological evolution for embodied cellular automata. Theoretical Computer Science, 2016; 633:19–27.

38.

Lipson

. Challenges and opportunities for design, simulation, and fabrication of soft robots. Soft Robot, 2014; 1(1):21–27.

39.

Sui

, Zheng

, Qi

, et al. Task-oriented hierarchical control of modular soft robots with external vision guidance. J Bionic Eng, 2022; 19(3):657–667.

40.

Peng

, Nabae

, Funabora

, et al. Controlling a peristaltic robot inspired by inchworms. Biomim Intell Robotics, 2024; 4(1):100146.

41.

, Awale

, Bacher

, et al. Scaling up soft robotics: A meter-scale, modular, and reconfigurable soft robotic system. Soft Robot, 2022; 9(2):324–336.

42.

Lee

, Her

, Jung

, et al. Snakeskin-inspired 3D printable soft robot composed of multi-modular vacuum-powered actuators. Actuators, 2023; 12(2):62.

43.

Zhou

, Xu

, Fang

. A CPG‐based versatile control framework for metameric earthworm-like robotic locomotion. Adv Sci, 2023; 10:2206336.

44.

Liu

, Zang

, Gao

, et al. Locomotion gait control of snake robots based on a novel unified CPG network model composed of Hopf oscillators. Rob Auton Syst, 2024; 179:104746.