Modeling and Fabrication of Soft Actuators Based on Fiber-Reinforced Elastomeric Enclosures

Unlike rigid actuators, soft actuators can easily adapt to complex environments. Understanding the relationship between the deformation of soft actuators and external factors such as pressure would enable rapid designs based on specific requirements, such as flexible, compliant endoscopes. An effective model is demonstrated that predicts the deformation of a soft actuator based on the virtual work principle and the geometrically exact Cosserat rod theory. The deformation process is analyzed for extension, bending, and twisting modules. A new manufacturing method is then introduced. Through any combination of modules, the soft actuator can have a greater workspace and more dexterity. The proposed model was verified for various fiber-reinforced elastomeric enclosures. There is good agreement between the model analysis and the experimental data, which indicates the effectiveness of the model.


Introduction
In the field of soft actuators, a fiber-reinforced elastomeric enclosure (FREE) developed by McKibben et al. [1][2][3][4][5][6] has attracted attention because of its simple structure, high forceto-weight ratio, and high pressure. These unique capabilities enable a variety of potential medical device applications [7,8]. In addition, soft actuators conform to the surroundings and move with dexterity [9][10][11]. Many soft actuators based on FREEs have been developed. They have a wide range of motions via combinations of programmed actions, including extension, bending, and twisting. Thus, it is essential to understand how to design a soft actuator to achieve a particular spatial motion. To make the design simpler and more systematic, there is a need to develop soft actuators that are easy to manufacture and program with respect to the complex environment.
As actuators based on FREEs can better match the behavior of biological muscles than others, research groups have begun to analyze the characteristics of soft actuators [12][13][14]. Finite-element analysis [15] is commonly used to optimize and guide the parameters for a soft actuator, and can also suggest designs. However, it can only analyze the movement of the foundation, and the calculations can be very large computationally.
Krishnan et al. [16,17] analyzed the mobility of generalized FREEs based on simple geometric relationships that result from the inextensibility of fibers and fluidic actuation. They examined the static behavior, such as extension, contraction, and twisting, by the principle of virtual work. However, this predicted only one-dimensional shape changes of soft actuators, without exploring actual application scenarios. To make it more flexible, some researchers arranged different fiber angles on the soft actuator to enable the single-drive actuator to perform a variety of actions [18][19][20]. These reports were based on geometric modeling, rather than a force analysis of the soft actuator.
Uppalapati et al. [21,22] combined multiple modules in parallel to maximize the workspace. They introduced an inverse formulation using Cosserat rods [23][24][25]. However, they did not establish a direct numerical relationship between pressure and deformation of a soft actuator. The actuator length was assumed to be constant, which is effective for bending or winding modules, but not for extension or twisting. Hence, there is a need to design a soft actuator that performs a specific motion, especially in complex environments.
To improve predictions, we used a geometric analysis model based on the principle of virtual work, combined with a mechanical analysis model based on Cosserat rod theory. Together, they predicted deformations of a soft actuator, especially gravity-induced shape deformation.
The advantages of using a combined model are: (1) The geometric analysis estimates the extension, contraction, and torsion angle at any time. The Cosserat rod theory cannot perform this requirement.
(2) The Cosserat rod theory analyzes soft actuator deformation via mechanical characteristics. External forces can be integrated into the model. The geometric analysis model cannot account for the effects of external forces.
Hence, the combination of these analytical models provides more accurate deformation features. Extension, bending, and twisting modules are fabricated, and deformations are verified. Experimental data agreed well with the proposed model. Finally, we randomly combined modules and performed a series of target capture experiments to demonstrate this approach.

Analytical Modeling of Actuator Segments
In this section, we describe in detail the Cosserat rod theory and the geometric modeling of soft actuators.

The Cosserat Rod Theory
The centerline curve of a rod with r(s) ∈ R 3 and material orientation R(s) ∈ SO(3), as functions of a reference parameter s ∈ [0, l], are the main characteristics, as shown in Figure 1.
drive actuator to perform a variety of actions [18][19][20]. These reports were based on geometric modeling, rather than a force analysis of the soft actuator.
Uppalapati et al. [21,22] combined multiple modules in parallel to maximize the workspace. They introduced an inverse formulation using Cosserat rods [23][24][25]. However, they did not establish a direct numerical relationship between pressure and deformation of a soft actuator. The actuator length was assumed to be constant, which is effective for bending or winding modules, but not for extension or twisting. Hence, there is a need to design a soft actuator that performs a specific motion, especially in complex environments.
To improve predictions, we used a geometric analysis model based on the principle of virtual work, combined with a mechanical analysis model based on Cosserat rod theory. Together, they predicted deformations of a soft actuator, especially gravity-induced shape deformation.
The advantages of using a combined model are: (1) The geometric analysis estimates the extension, contraction, and torsion angle at any time. The Cosserat rod theory cannot perform this requirement.
(2) The Cosserat rod theory analyzes soft actuator deformation via mechanical characteristics. External forces can be integrated into the model. The geometric analysis model cannot account for the effects of external forces.
Hence, the combination of these analytical models provides more accurate deformation features. Extension, bending, and twisting modules are fabricated, and deformations are verified. Experimental data agreed well with the proposed model. Finally, we randomly combined modules and performed a series of target capture experiments to demonstrate this approach.

Analytical Modeling of Actuator Segments
In this section, we describe in detail the Cosserat rod theory and the geometric modeling of soft actuators.

The Cosserat Rod Theory
The centerline curve of a rod with   The v(s) and u(s) are kinematic variables that represent the linear and angular rates of change of r(s) and R(s) in a local frame, respectively. Thus, they have the following relationship: . .
Actuators 2021, 10, 127 3 of 12 whereû(s) represents the skew symmetric matrix of the vector u(s). This is defined: From the force and moment balance analysis of the rod in Figure 1, the following can be obtained: where n(s) indicates the internal force and m(s) indicates moment vectors. R(s)v(s) indicates the tangential direction of a point on the soft actuator. f d (s) is the distributed force acting on the actuator and L(s) is the applied moment distribution. We use linear constitutive laws to map these variables to the internal forces and moments, as shown in Equations (5) and (6): where is the stiffness matrix for shear and extension, and K bt (s) is the stiffness matrix for bending and torsion. G is the shear modulus, and I xx , I yy , and I zz are the second moments of area for the tube cross sections about the principal axes.
It is necessary to make a brief explanation here. For the bending and twisting modules v(s) = v 0 (s) = 0 0 1 T and the axial force F e = 0. In other cases, v(s) = v 0 (s) = 0 0 λ 1 T .
The explicit model equations can be obtained by combining Equations (1)-(6) and the following: .
Then, the boundary conditions are determined. At Matlab ode45 was used with the shooting method to solve the boundary-value problem.

Geometric Model Analysis
The relationship between soft actuator deformation and the fiber angle is given by: where θ = l 0 tan(α) , λ 2 = r * r 0 , and λ 1 , λ 2 indicate axial elongation and radial elongation, respectively, and α, β ∈ 0, ± π 2 indicate the helical angles of the fiber, as shown in Figure 2. δ is the net angular deformation of the actuator's end. l 0 and r 0 are the initial actuator lengths. l * , r * are the deformed elemental axial length and radius, respectively. , and 1 λ , 2 λ indicate axial elongation and radial elongation, respectively, and α , indicate the helical angles of the fiber, as shown in Figure 2. δ is the net angular deformation of the actuator's end. 0 l and 0 r are the initial actuator lengths. * l , * r are the deformed elemental axial length and radius, respectively. The fibers were inextensible, so the relationship between input pressure and axial elongation can be obtained by the virtual work principle, as follows: The fibers were inextensible, so the relationship between input pressure and axial elongation can be obtained by the virtual work principle, as follows: where P is the pressure and V is the volume of the actuator cavity. V t is the volume of rubber, and E is Young's modulus of the material. This was numerically inverted to solve for λ 1 , given P from Equations (9)-(13).

Modeling of Extension, Bending and Twisting
To model extension and bending, the helical angles of the fibers were set to α = −β, (α, |β|) > 54.73 • . If the angles are symmetrical, the soft actuator becomes longer under continuous air pressure and will not twist. Therefore, Equations (9)-(10) can be simplified as follows: The only difference between the bending and extension modules is that the bending module must add another fiber as the binding layer. The added fiber allows the soft actuator to bend towards the fixed side. The process of change with the bending module is shown in Figure 3.
The twisting module is a combination of the bending and rotating modules, as shown in Figure 4. The helical angle of the fiber had an asymmetrical arrangement, with a motion similar to the predator action of an octopus arm. In this case, Equation (14) is not applicable and λ 1 can only be solved with Equations (9)- (13). Finally, the relationship between pressure and the deformation of the soft actuator can be obtained through the principle of virtual work under ideal conditions. All of the above cases were obtained for an ideal state, without considering the effects of external forces, such as gravity. In Figure 5, the gray cylinder on the left depicts the length of an actuator without the effect of gravity. The dark gray cylinder on the right depicts the change in length when considering the effect of gravity. Hence, there is a large difference when gravity is considered for a soft actuator. The only difference between the bending and extension modules is that the bendin module must add another fiber as the binding layer. The added fiber allows the soft actu ator to bend towards the fixed side. The process of change with the bending module shown in Figure 3. ( ) The twisting module is a combination of the bending and rotating modules, as show in Figure 4. The helical angle of the fiber had an asymmetrical arrangement, with a motio similar to the predator action of an octopus arm. In this case, Equation (14) is not applic ble and 1 λ can only be solved with Equations (9)- (13). Finally, the relationship betwee pressure and the deformation of the soft actuator can be obtained through the principle virtual work under ideal conditions.   The twisting module is a combination of the bending and rotating modules, as show in Figure 4. The helical angle of the fiber had an asymmetrical arrangement, with a motio similar to the predator action of an octopus arm. In this case, Equation (14) is not applica ble and 1 λ can only be solved with Equations (9)   All of the above cases were obtained for an ideal state, without considering the effects of external forces, such as gravity. In Figure 5, the gray cylinder on the left depicts the length of an actuator without the effect of gravity. The dark gray cylinder on the right depicts the change in length when considering the effect of gravity. Hence, there is a large difference when gravity is considered for a soft actuator. To describe the above actions more accurately, the Cosserat rod theory was used to improve the soft actuator model under external loads. Initially, the equivalent axial force generated by internal pressurization of the tube was determined with Equation (16). This equation makes a preliminary hypothesis for the model derivation.
where t A is the cross-sectional area of the soft actuator and is equal to To describe the above actions more accurately, the Cosserat rod theory was used to improve the soft actuator model under external loads. Initially, the equivalent axial force generated by internal pressurization of the tube was determined with Equation (16). This equation makes a preliminary hypothesis for the model derivation.
where A t is the cross-sectional area of the soft actuator and is equal to A t = π r 2 * − r 2 0 .

Manufacturing Soft Actuators
The actuators needed five materials: a hollow latex tube, fibers, rubber cement, a latex tube surface-treatment agent, and a metal rod, as shown in Figure 6a,b. The latex tube had 15-mm inner and 20-mm outer diameters. To describe the above actions more accurately, the Cosserat rod theory was used to improve the soft actuator model under external loads. Initially, the equivalent axial force generated by internal pressurization of the tube was determined with Equation (16). This equation makes a preliminary hypothesis for the model derivation.
where t A is the cross-sectional area of the soft actuator and is equal to

Manufacturing Soft Actuators
The actuators needed five materials: a hollow latex tube, fibers, rubber cement, a latex tube surface-treatment agent, and a metal rod, as shown in Figure 6a,b. The latex tube had 15-mm inner and 20-mm outer diameters.  The soft actuators were fabricated by mounting the latex tube on the metal rod, which had a diameter equal to the inner diameter of the latex tube. To paste the fiber and latex tube together, a layer of treatment agent (J-740, DOLDSUN TIME) was applied to the surface of the latex tube. This agent enhanced the fiber restraint effect.
Then, the metal rod was fixed for machine control of the fiber-winding angle on the latex tube. The molding effects are shown in Figure 7. Finally, we used rubber cement to adhere the fiber to the latex tube. The soft actuators were fabricated by mounting the latex tube on the metal rod, which had a diameter equal to the inner diameter of the latex tube. To paste the fiber and latex tube together, a layer of treatment agent (J-740, DOLDSUN TIME) was applied to the surface of the latex tube. This agent enhanced the fiber restraint effect.
Then, the metal rod was fixed for machine control of the fiber-winding angle on the latex tube. The molding effects are shown in Figure 7. Finally, we used rubber cement to adhere the fiber to the latex tube.

Experimental Validation
To fabricate extending, bending, and twisting actuators, we used a hollow latex tube with  Table 1.

Experimental Validation
To fabricate extending, bending, and twisting actuators, we used a hollow latex tube with E = 7 × 10 5 N/m 2 and G = 2.33 × 10 5 N/m 2 . The fiber reinforcement was Kevlar. It was assumed that the fiber was inextensible. For the extension and bending modules, the fiber angles were set to α = −β = 85 • , and α = 85 • , β = −80 • for the twisting module. Detailed parameters are listed in Table 1.

Module Testing
The modules were verified separately with experimental data. The model data of the combined virtual work principle and the Cosserat rod theory were compared and analyzed with experimental data.

Extension Module Test Experiment
The length of the extension module was set to 0.478 m. In Figure 8, the red line is the experimental data for the actuator length, and the blue line is the length of the extension module predicted by the model. The experimental data were always greater than the data predicted by the model, and the maximum error was 0.0056 m. The reason for the large error may be because the extension module was not an ideal straight line during the fabrication process. There was a certain degree of bending. Therefore, with continuous pressure, the soft actuator was slightly bent along the elongation, which led to the actual length of the soft actuator always being less than the length predicted by the model. In general, the data predicted by the model matched the experimental results well, and were effective enough to predict the shape change of the model.

Bending Module Test Experiment
The length of the bending module was 0.378 m. Then, with input pressures P from 40 kPa to 110 kPa, we used the analytical model to predict the shape of the actuator, as shown in Figure 9. The direction of gravity was positive along the z-axis. The abscissa represents the deviation of the soft actuator from the y-axis direction. Experimental data under corresponding conditions were also acquired. The position of the soft actuator must

Bending Module Test Experiment
The length of the bending module was 0.378 m. Then, with input pressures P from 40 kPa to 110 kPa, we used the analytical model to predict the shape of the actuator, as shown in Figure 9. The direction of gravity was positive along the z-axis. The abscissa represents the deviation of the soft actuator from the y-axis direction. Experimental data under corresponding conditions were also acquired. The position of the soft actuator must be determined via image processing. The possible reason for the error is that the curvature of the soft actuator was not a constant curvature model with increasing input air pressure. Thus, there may have been deviations in the estimated position of the actuator. Another reason is that the cover at the end of the actuator had weight, which could make the predicted shape deviate from the actual shape. In Figure 9, there is good agreement between the analytical model and the experimental data. To better evaluate the accuracy of the analytical model, we used end-position errors to measure the model. Thus, the prediction model was evaluated by determining the position error of the end point between the model and experiment. The error was calculated by: The possible reason for the error is that the curvature of the soft actuator was not a constant curvature model with increasing input air pressure. Thus, there may have been deviations in the estimated position of the actuator. Another reason is that the cover at the end of the actuator had weight, which could make the predicted shape deviate from the actual shape.

Twisting Module Test Experiment
The length of the twisting module was 0.7 m. As its motion was three-dimensional, a helical radius was used to measure the model, extract the experimental data, and evaluate the accuracy of the model. The radius could be obtained by the virtual work principle and the Cosserat rod theory. Initially, the λ 1 of the actuator under ideal conditions was solved by using the virtual work principle. Then, the Cosserat rod theory was used to obtain the actuator attitude information. The curvature of the actuator centerline could be obtained by this method. The relationship between curvature and the helical radius could be obtained from previous reports [22]: where τ is the twist per unit length. ξ represents the angle between the tangent direction at the end of the actuator and the z-axis direction of the global coordinate system.
The change in radius of the twisting module could be obtained by the above method, as shown in Figure 10. The experimental data were always larger than the predicted data from the model. The maximum error was 0.01927 m, which may be because the cover at the end of the torsion module was made of metal and had a certain weight that affected the helical radius. However, the agreement between the spatial twisting model and experimental data was quite high, indicating the reliability of the proposed model.

Model Series Motion Experiment
The function of a single module is limited. Therefore, to highlight the operation of the soft actuator, we fabricated the extension, bending, and twisting modules at different scales. The twisting module was used as the end actuator, and it captured the target, as shown in Figure 11a. Then, we combined the torsion module with the bending and elongation modules to grasp the target. In Figure 11b,c, the twisting, bending, or extension modules could grasp the target. as shown in Figure 10. The experimental data were always larger than the predicted data from the model. The maximum error was 0.01927 m, which may be because the cover at the end of the torsion module was made of metal and had a certain weight that affected the helical radius. However, the agreement between the spatial twisting model and experimental data was quite high, indicating the reliability of the proposed model.

Model Series Motion Experiment
The function of a single module is limited. Therefore, to highlight the operation of the soft actuator, we fabricated the extension, bending, and twisting modules at different scales. The twisting module was used as the end actuator, and it captured the target, as shown in Figure 11a. Then, we combined the torsion module with the bending and elongation modules to grasp the target. In Figure 11b,c, the twisting, bending, or extension modules could grasp the target.

Conclusions
Although there have been many reports on model analyses for soft actuators, especially for FREEs, few systematically analyzed separate modules, especially considering the effects of external forces. Here, we improved the description of spatial state changes in soft actuators. Firstly, a geometrically exact steady-state model (Cosserat rod theory) and the virtual work principle were used to predict deformation characteristics. The most important contribution is that the combination of these two models increased the design speed of the soft arm and improved the accuracy of shape prediction. Secondly, to verify the effectiveness of the model, we used the model to predict extension, bending, and twist-

Conclusions
Although there have been many reports on model analyses for soft actuators, especially for FREEs, few systematically analyzed separate modules, especially considering the effects of external forces. Here, we improved the description of spatial state changes in soft actuators. Firstly, a geometrically exact steady-state model (Cosserat rod theory) and the virtual work principle were used to predict deformation characteristics. The most important contribution is that the combination of these two models increased the design speed of the soft arm and improved the accuracy of shape prediction. Secondly, to verify the effectiveness of the model, we used the model to predict extension, bending, and twisting deformation actuators. Then, a comparison of the model analysis and the experimental data for the three modules verified that the result errors were 0.0056 m, 0.0519 m, and 0.01927 m. Finally, to expand the function of the soft actuator, the extension, bending, and twisting modules were randomly combined to grasp objects in space by controlling the pressure of the three modules.
A soft actuator constructed using FREEs can easily change its shape by means of external forces. Hence, in future work, we plan to change this structure into a module with variable stiffness, then we will optimize the structure of the soft actuator to minimize the volume. Finally, we plan to examine dynamic modeling analysis of the soft actuator because most modeling of the flexible actuator is still based on quasi-static motion.
Author Contributions: X.H. was the supervisor providing funding and administrating the project, and he reviewed and edited the manuscript. Z.C. carried out the investigation, methodology, analysis, and validation, and wrote the original manuscript. A.Z. and T.L. have made great contribution to the application of MATLAB in simulation. Z.Q. built a platform for the verification model and helped to make different soft actuator modules. S.L. helped to test the parameters of the soft actuator, and also helped to perform some verification experiments. All authors have read and agreed to the published version of the manuscript.