Finite Element Method Analysis and Structure Design of Stiffness-tunable Beam-shaped Material

Most surgical tools are made of metal or soft rubber, the stiffnesses of which may not be suitable for some surgeries. There is a need for stiffness-tunable materials for use in surgical tools or a mechanism for tuning the stiffness of a material without a chemical or thermal reaction. We have proposed a stiffness-tunable beam-shaped material with pneumatic access for stiffness control. We have also shown its application in grasping a soft organ stably during laparoscopic liver surgery. Its stiffness became 5.47 times greater by inducing -80 kPa negative pressure inside the device. The beam structure had a wavy composition of hard and soft rubbers. In this paper, we evaluate the finite element method (FEM) performance to determine the optimal design of the stiffness-tunable beam. FEM resulted in a 11.3-fold improvement of stiffness tuning for wavy shapes, that is, sequentially aligned pillars with an adequate pitch, and increased the stiffness, one of the shaping parameters, 11.8-fold.


Introduction
Grasping and stably holding the liver is a significant issue in surgical safety as the surgeon's handling force and torque must be precisely transmitted to the organ in laparoscopic liver surgeries. Inappropriate manipulations can harm liver tissues or blood vessels, resulting in a rough resection line or serious bleeding. In present clinical procedures, a surgeon grasps and pulls the liver using a grasper at a single contact point, aligns the liver, and then applies adequate tension to maintain a smooth transection. This procedure may cause excessive grasping pressure to be applied to the liver and hence tear the tissue. To avoid such tissue damage, some surgeons employ retractors to form a surgical contact area in resection. However, slippage or excessive grasping pressure on the liver may occur. In addition, to avoid the risk of bleeding, care must be taken that the retractors do not come into contact with the cut face of the liver. Suction fixation may be a solution to these problems. (1)(2)(3)(4) Multiple suction pads adhering to the liver surface can achieve higher friction resistance. Their adequate spatial distribution can also reduce the contact pressure at each pad. Since manufacturing a patient-specific shaped manipulator with a rigid structure is impractical, a manipulator that can be instantaneously adapted to an arbitrary organ surface while maintaining high rigidity would be ideal. Such a manipulator requires both shape followability and position-holding ability, which may be in a tradeoff relationship with manipulator stiffness. One promising way to overcome this problem is to adopt a stiffnesstuning mechanism. This mechanism first reduces the stiffness to enable the manipulator to follow the liver surface, and then increases it to hold the liver in position. Some researchers have used wire tension in rigid links for stiffness tuning. (1)(2)(3) However, wire tension resulted in friction between rigid links and prevented their rotation. Thus, the whole link structure was in a high-stiffness state. On the other hand, miniaturization is challenging because high wire tension and structural complexity require the link mechanism to be of a certain size. We previously developed a stiffness-tunable beam-shaped material and applied it to assist liver manipulation in laparoscopic liver surgery, as shown in Fig. 1. (5,6) The stiffness of the device was increased 5.47fold through automatic pneumatic control of the device interior. In addition, we showed the stiffness-tuning principle and the effectiveness of the wavy-bounded design of a hard and soft material composite. However, we did not address any strategy of device design or parameter tuning. In this study, we employ the finite element method (FEM) to compute the stress distribution of objects and use the results to design an optimal device structure. FEM analysis is commonly used to optimize object structures. It has often been employed in designing buildings, machines, and other infrastructures. In addition, in recent decades, it has been used to optimize material structures, which is the essence of the research on material informatics. (6)(7)(8)

Device principle and configuration
The principle is given to describe simply how the device transforms to a hard stiffness state when using a homogeneous beam material. Figures 2(a) and 2(b) show the cross-sectional profiles of the device. The device consists of two soft beams, a chamber sponge, a membrane  covering the device, a supplementary extrusion wing membrane, a vacuuming tube, and grasping bars. The two soft beams of the device are nonstretchable but easy to bend. They are made of rubber that contains nonstretchable fibers. The device is connected to a medical vacuuming pump to remove the underlying external air. As shown in Fig. 2(a), the two soft beams are free to slide in this state, which is the low-stiffness state. After the device adhered to the liver, as shown in Fig. 2(b), the device bottom was covered completely by the liver and the inner air was exhausted, leading to a transition of the device stiffness owing to fixation of the two soft beams. Let b and h be the width and height of the cross section of a single beam, respectively, as shown in Fig. 2(c). The bending stiffness of a square beam is where E is the elastic modulus and I single_beam is the second moment of area for one beam. Before the air is exhausted, the two beams are free to slide against each other. Thus, the device stiffness is double the single-beam stiffness and can be described as (2) After the air inside the device is exhausted, the two beams are pressed together by atmospheric pressure from the outside, and the beams can be considered to be unified into a double-thick beam. Its height becomes 2h, as shown in Fig. 2(d), and thus the stiffness is (3) Because EI constrained_beams is fourfold greater than EI unconstrained_beams , the pneumatic operation can change the device stiffness. Details of the beam design are described in the next section. The sponge chamber at the bottom stabilizes the airway for suction of the liver surface. It is covered with soft gauze on the bottom side to reduce liver damage caused by suction and pressure. The extrusion wing around the bottom edge is a thin membrane made of silicone and adheres to the liver surface and acts as a seal. A tube to exhaust air is attached at the top of the device.

Theoretical analysis of stress distribution
The stress distribution of the beam is given following the condition that the device adheres to the liver surface and is pulled up at a point on the beam. Figure 3(a) shows the dimensions of the square beam, the two constrained beams, and its coordinate system. The directions of the x-and y-axes are horizontal and vertical, respectively. The y-axis is positive downward. The length, width, and height are l, b, and 2h, respectively. The device was pulled up at x = 0 with force F. The external force applied by the liver to the device was uniformly distributed horizontally at the bottom. Let uniformly distributed load be w and the left-and right-side lengths be l 1 and l 2 , respectively. Torque, N = w(l 2 2 − l 1 2 )/2, given at x = 0 cancels the device rotation caused by the unbalance of l 1 and l 2 . The beam material is assumed to be homogeneous. Under this condition, bending stress is Shear stress is derived as Equation (5) gives a parabolic distribution for y with absolutely maximized shear stress at (x, y) = (0, 0), as shown in Fig. 3(c). The bending stress distribution, expressed by Eq. (4), means that the tolerance for tension and compression should be ensured at y = ±h., i.e., at the top and bottom faces of the beam. In addition, the shear stress distribution, expressed by Eq. (5), means that the device should be resistant to shear stress at y = 0, i.e., at the boundary face between the two pressed beams. The device suction performance at the liver surface is pbl, which is determined by the exhaustive pressure p and the bottom suction area of the device, bl.
Considering these conditions, beams with one wavy side, as shown in Fig. 4, were employed to make the stiffness-tuning mechanism effective. The beam bases, which are colored black, were made of the resin 60S3M210 (BANDO Chemical Industries, Ltd., Japan) and the gaps were filled with silicone. Silicone sheets (Togawa Rubber Inc., Japan) were affixed to the internal side of each beam to prevent slipping between beams. The wavy shape of the beams worked well to allow easy bending and anchored the silicone sheets to prevent plastic deformation at the medial boundary between the two beams. The distortion resistance of beams under the bending condition was due to the upper beam being curled while the lower beam was outwardly bent. The beam bases were subject to major tension and compression at the top and bottom, respectively. The convexity of the base enabled it to resist the sliding force from the medial rubber sheets. Soft silicone filling the base gaps allowed beam bending and prevented wavy deformation of the medial rubber sheets. Five vertical holes were formed in each beam to ensure an airway. In addition, a runnel was placed at the center along the long axis of the beam. A handlebar was attached to the bottom beam to support easy grasping with laparoscopic forceps.

FEM stress computation and optimization of beam structure
The stress distribution of the stiffness-tunable beam was computed by FEM. The tested beam models were of three types: homogeneous, with hill-shaped pillars, and with exponential-shaped pillars, as shown in Fig. 5. The hill-shaped pillar had the same shape as a timing belt driving a gear wheel. The exponential-shaped pillars might be considered upright but all crutch points were made thicker. The device width, height, and length were 8.0, 4.88, and 50.0 mm, respectively. Young's moduli, E 1 and E 2 , were 20.0 and 3.0 MPa for high-and low-rigidity materials, respectively. Their Poisson's ratio was 0.49. High-friction rubber sheets were placed between these beams. Their Young's modulus and Poisson's ratio were 210 × 10 3 MPa and 0.3, respectively. Mesh elements used in FEM were 0.5-mm-size quadratic tetrahedrons. The numbers of nodes and elements were respectively 249480 and 171651 for homogeneous types, 1721301 and 1040611 for hill-shaped-pillar types, and 1150355 and 710832 for exponentialshaped-pillar types. A vertical force of 0.1 N was loaded to the models in each test. Fusion 360 CAD (Autodesk Inc., USA) was used for FEM computation. To compare the three beam shapes, stress distributions were evaluated by the three-point bending test as shown in Fig. 6(b). The stress profiles at some cross-sectional lines were given in this test. Then, the stiffness transition was evaluated with cantilever bending simulation to optimize the beam design, as shown in Figs. 6(c) and 6(d). The model with exponential-shaped pillars was used in this test. As shown in Figs. 6(e) and 6(f), pillar pitch was varied as 1.5, 3.0, 6.0, and 9.0 mm. The pillar occupancy was changed from 10 to 90 in 10% increments. Beam-body thickness, which is the inverse of gap depth occupancy, was changed from 10 to 90% in 10% increments. Let F, l, and δ be load, longaxis length of the beam, and vertical deformation displacement, respectively. Bending rigidity is then given by  As well, bending rigidity ratios are given with bending rigidities in low-and high-rigidity states as where δ 1 and δ 2 are vertical deformation displacements of the upper and lower beams of the device.

Results
Stress distribution of the beams was computed by FEM as shown in Fig. 7. Figures 7(a), 7(c), and 7(e) show normal stress, which is tensile or compressive stress, and Figs. 7(b), 7(d), and 7(f) show shear stress. Figures 7(a) and 7(b) are for homogeneous-material beams. The stress distributions correspond to the theoretical results shown in Figs. 3(b) and 3(c), respectively. In the horizontal direction, stress is maximum at the loading point. In the vertical direction, normal stress is maximum at the top and minimum at the bottom. Shear stress is maximum at the vertical middle.   shown in Fig. 9, the beams with hill-shaped and exponential-shaped pillars gave the same value of 11.3 for bending rigidity gain. For homogeneous-material beams, bending rigidity gain was 3.95, which is similar to the theoretical value. Low and high bending rigidities and their gains are given with respect to each pillar pitch a, pillar occupancy α at each pillar pitch, and beam body occupancy β in the entire thickness, as shown in Fig. 10. Bending rigidity did not show any obvious trend but gain showed a slight decline with pillar pitch. Rigidity gain showed an 11.4-fold increase between pillar occupancy α of 40 and 50% in Fig. 10(d). This was due to the increase in bending rigidity in the low-rigidity state (triangles), as shown in Fig. 10(c). Between 40 and 50%, the normal stress distribution changed broadly. One of the reasons behind this might be that the bending rigidities of the beam body and soft sheet attached to the beam bottom became equal. Figures 11(a) and 11(b) show details of the stress distribution change. When α exceeded 50%, the beam body was harder than the soft sheet, causing the stress to be distributed in the beam body; thus, the stress in the soft sheet could be ignored. When α was less than 40%, rigidities of the beam body and soft sheet became equal. This caused the singularity shown in Fig. 10(d). In Figs. 10(e) and 10(f), bending rigidity showed the peak of 11.8-fold increase between beam-body thickness occupancy β of 40 and 50%. Figures 11(c) and 11(d) show detailed stress profiles of the downward inflection. At 40% pillar occupancy, the normal stress changed from the top and bottom to the entire beam as a result of the equalization of the bending rigidities of the beam body and the soft sheet attached to the beam bottom. At 50% pillar occupancy, the stress showed a major change between the top and bottom of the beam body because the beam-body rigidity was high enough for normal stress to be concentrated in the beam body.

Discussion
FEM analysis enabled us to quantify and visualize the stress distribution of stiffness-tunable beam materials, as shown in Fig. 7. It showed that a wavy structure of rigid and soft rubbers could enable effective tuning of the stress distribution. As shown in Fig. 8, FEM analysis showed that the maximum stress appeared in the rigid material. The wavy structure was effective for optimal stiffness tuning.
The beam having exponential-shaped pillars showed the best performance with an 11.3-fold increase in stiffness, as shown in Fig. 9. Exponential but almost-straight pillars are often seen in designs of buildings and bridges. Concrete blocks might also be similarly shaped. It may be one of the optimal shapes in material mechanics. Compared with the rigidity gain of 5.48 reported in our previous paper, (6) FEM analysis gave a rigidity gain of 11.3, which is approximately twice as high. The reason might be that irregular stress unbalance and stiffness-sheet distortion, which occur in reality, such as sheet wrinkles, were not considered in the FEM simulation. Although the absolute value of beam stiffness in the simulation is higher than the true value, and we should also compare simulation results with real measurements in the future, the trend of beam rigidity change was captured in beam shape designing with FEM simulation. In addition, we assumed that the beam materials were homogeneous but did not have any heterogeneous stretching characteristics. This assumption is a possible cause of some of the differences between real and simulation results.  Figures 10(d) and 10(f) showed that the inflection points of curves were in the midrange of pillar occupancy rate α and beam body occupancy rate β. These were caused by the decrease in high-stiffness rigidity and the increase in low-stiffness rigidity with respect to α or β, as shown Figs. 10(c) and 10(d). In these results, we observed points of singularities in the change of rigidity gain with respect to α and β. As shown in Fig. 11, they appeared when the rigidities of the beam body and soft sheet became equal. The center of second-area moment, y 0 , can be approximately given by where E 1 and E 2 are Young's moduli for high-and low-rigidity materials, respectively. They were 20 and 3 MPa, respectively, for the materials we used. Young's modulus of the silicone filling the beam gaps was ignored in this numerical analysis. y 0 was plotted against α with β = 40% and against β with α = 40%. In our approximate theoretical analysis, the singularities occur at 49% for α and at 38% for β, as shown in Fig. 12. This means that the center of the second-area moment shifts into the area of the beam body from below. Although the singularity point depends on the rigidity moduli of the beam body and soft sheet, the observation of this phenomenon might be meaningful in designing stiffness-tunable materials.
In this paper, we designed beam shapes for a stiffness-tunable device to apply it to liver manipulation in laparoscopic surgeries with optimized beam shape parameters α and β. More precise optimization of shape parameters will enable complex device geometries, which would then require the introduction of additional manufacturing equipment such as 3D printers to shape such devices.

Conclusions
FEM analysis reproduced the trend of device stiffness change, indicating the feasibility of its use in device shape optimization. A wavy shape of the beam of surgical tools resulted in optimal device structure parameters that enabled stiffness tuning at both low-and high-rigidity levels. The optimal shape of the beam had pillar occupancy and beam body occupancy rates at the midrange of the range of values obtained in the results.