Figure 2 shows the schematic diagram of the experimental setup. In Figure 2A, the motor block drove the flexible shaft (a super-elastic nitinol tube with OD/ID: 1.27/0.97mm, NI207130, GOODFELLOW®) into the curved channel (glass fiber reinforced plastic tube OD/ID: 3/1.5 mm). The motor (2342S024C, FAULHABER®) and the encoder (2RM3600-D, SCANCON®, 3,600 pulse/rev) clamped the flexible shaft by two rubber rollers of 13.5 mm diameter. A housing tube was used to support the remaining nitinol tube outside the motor block. It also reduced the vibration of the nitinol tube. The overview of the testbed is shown in Figure 2B. The glass fiber reinforced plastic tube (channel) was bent into desired shapes and fixed on to the perforated board. The top-view image of the curved channel with the checkerboard was taken. The checkerboard was used as a reference when calculating the curvature along the channel. In Figure 2C, the motor block was mounted on a linear slider (SSEBWZ16-110, MISUMI®). Two loadcells (FUTEK® LSB200, capacity 4.5 N, rated output 0.1 mV/V) sandwiched the motor block were fixed on two linear sliders. Tuning these linear sliders could adjust the pretension applied on the motor block.

MATLAB SIMULINK interfaced with a QPIDe data acquisition board from Quanser Quarc® was used to control the motor and acquire and synchronize sensor readings. The motor was controlled to insert the flexible shaft into the channel at a speed of 20.0 mm/s. Once insertion length s = 700 mm was reached, the motor would rest for 5 seconds and then withdraw the shaft to the starting point. The motor rested for 5 seconds again and this process was repeated two times.

A force diagram is provided in Figure 2D. The insertion force required is denoted as   F i n s   . Both load cells are pre-loaded to   F 0   before measurement. The load cell readings are denoted as F 1 and F 2 . The forces from load cell 1 to the motor block is   Δ F 1 = F 1 − F 0 , and the force from load cell 2 is   Δ F 2 = F 2 − F 0 . The total insertion force is as follows:   F i n s =   Δ F 1 −   Δ F 2 = F 1 − F 2

A second-order low-pass filter (cutoff frequency, sampling rate) was applied for both load cells. Data recorded from insertion distances 30–600 mm was used for this study. Then, the insertion length L   =   570  mm and s ∈ [ 0,570 ]  mm. In this range, the motor block was in its most stable condition and the shaft is driven with constant speed. After calibration of all load cells, the platform was used to collect data. Figure 3 shows three repeated measurements on a channel, and the maximum deviation is 3% comparing with the mean value. It indicates the insertion force for a given shape is consistent, and the measurement system is reliable. Note that the frictions of the housing tube and linear bearing are neglectable.

The top-view image of the channel and the checkerboard was taken. The point clouds of the channel could be obtained via imaging processing in terms of the pixels. These pixels could be converted to real-world coordinates with respect to the checkerboard. Then, a parameter curve could be fitted to the point clouds and then used to calculate coordinates, curvatures, and accumulated bending angles.

First, the camera parameters were obtained using the checkerboard fixed on the platform using MATLAB camera calibration function “estimateCameraParameters”. Then, the object in the image could be measured using the checkerboard and camera parameters. To illustrate the accuracy of this measurement, two checkerboards were placed diagonally on one A4 paper, as shown in Figure 4A. One checkerboard was the reference checkerboard, and the other was the targeting/object checkerboard. Figure 4B shows the estimations and ground truth of the objects. The average and maximum distances between the estimates and ground truth are 1.26 and 1.88 mm, which are small enough for this application.

Then, we took the image of the channel during the insertion experiment. These images were post-processed to obtain the channel shape. Figure 5 shows the steps of extracting the channel coordinates with respect to the checkerboard from images taken. First, the fisheye effect (lens distortion) was removed using the camera parameters. Then, as shown in Figures 5A, a greyscale image was obtained by saturating the top and bottom 1% values of the green channel of the image. Then the greyscale image was converted to a binary image with a threshold of 180 as Figure 5B. Then the binary image would be converted to point clouds in coordinates form. The noises were removed from the tube objects, resulting in a clean background in Figure 5C. The point clouds were fitted with 1) straight lines for the beginning and ending 2) ninth-order polynomial for the curved section. A fifth-order spline interpolation replaced the intersection zone of linear fittings and the polynomial fitting to obtain a smooth curvature profile as shown in Figure 5C. The starting point of the insertion shaft is set to be the origin. Figure 5D shows the point cloud obtained and the fitting curve. We obtained the acutal shape P c h a n n e l = [ x , y ] by converting these coordinates in the pixel to mm using the checkerboard. The origin of the acutal shape P c h a n n e l was fixed as shown in Figure 5A. All negative x points of the curve were discarded. The insertion length of m^th point in the curve is calculated as the length of the fitting curve from the original point to the m^th point: s m = {                                                                                     0                                                                                               f o r   m = 1 ∑ i = 2 m ( x i − x i − 1 ) 2 + ( y i − y i − 1   ) 2                             f o r   m > 1

Then curvature was calculated based on the fitting curve using LineCurvature2D (Kroon 2021). The incremental angle for each point is calculated by d θ i = {                                                                   0                                                           f o r   i = 1   o r   i = n tan − 1   Δ y i + 1 Δ x i + 1 − tan − 1   Δ y i Δ x i                   f o r i ≠ 1   a n d   i ≠ n w h e r e   Δ x i =   x i − x i − 1   a n d   Δ y i =   y i − y i − 1    

The accumulated bending angle from beginning to m^th point is calculated by: θ m =   ∑ i = 0 n   | | d θ i | |   .

The insertion forces were measured when inserting the shaft into the curved channels. The channel curvatures were measured from the top-view image. The insertion force and respective channel curvatures were grouped as one set of data. By changing the shape of the channel, different insertion forces were obtained. As summarized in Table 1, we collected the channel curvatures and insertion forces of 120 one-bend shapes (indexed from #1 to #120) and 39 two-bend shapes (indexed from #121 to #159). One-bend shapes #35, #41& #106 in Figure 6A and two-bend shapes #141, #145 & #159 in Figure 6B are shown as an example The ranges of curvature and accumulated angle are also listed in Table 1. Moreover, the curvature was set to be positive for one-bend shapes. For two-bend shapes, the first bend had positive curvature, but the second bend could be positive or negative. The insertion force profiles and curvatures for one-bend shapes are shown in Figures 6B,C and those for two-bend shapes are shown in Figures 6E,F. It can be observed that the insertion force changes in a nonlinear manner in Figure 6C as the curvature changes in Figure 6B. Then, these shapes were grouped into two different datasets for the neural network model as Table 2. We used 105 one-bend shapes as training/dev set and 15 one-bend shapes as test set. Dataset 2 mixed one-bend shape and two-bend shapes. All one-bend shapes and 32 two-bend shapes were in training/dev set while 7 two-bend shapes were in test set.

Then, we used the data driven method to investigate association between the insertion force profile and the channel shape.

First, the data obtained (insertion force profile and channel shape) are sequential. At any time, the curvature along the channel could be treated as function of insert length, provided the clearance is small. Then we also record the insertion force profile with respect to insertion length. The sequential data could be handled by the long-short-term-memory (LSTM) network. The bi-directional LSTM is formed by stacking two LSTM of opposite directions. The bi-directional LSTM has advantages in processing the context of a sequential data.

The structure of the proposed deep neural network model is shown in Figure 7. The network has an input as: F i n s ( s i ) ( k )     at every insertion length step i for shape k. The label or output is the curvatures   κ ( s i ) ( k )   at insertion length step i for shape k. The input data was first fed to the bi-directional LSTM layers and then output through a MLP network. A simple way of enhancing the capability of a LSTM network is to increase the number of layers. However, if only LSTM is used, as the number of layers increases, more variables inside the LSTM need to be handled, and the training time will significantly increase. Therefore, multiple bidirectional LSTM layers are used for feature extraction and multilayer perceptron (MLP) is stacked to handle the information from bidirectional LSTM, which has been proved its advantages in (Zhao et al., 2017).

The network parameters were eventually selected using Sklearn. gridsearchCV: four layers of bidirectional LSTM, two layers of MLP with dropout layers in between, 120 hidden neurons for all layers. The dropout ratio was fixed at 0.5. Adam as optimizer and Huber loss-function are the best pick of parameters for the proposed neural network. The training stops at 5,000 epochs, where loss has converged. The configuration of the workstation: processor: Inter® Core (TM) i5-7500; RAM: 16G; GPU: Nvidia GTX 1050Ti. The CuDNNLSTM (Keras.layers) was used instead of normal LSTM layers to reduce training times. With parameters listed above, 5-fold cross validation was performed, i.e., the proposed neural network was trained and tested on five splits of Dataset 1. The training scores and testing scores for each split dataset were calculated. Finally, it is to check whether the network fits all the split datasets equally.

Table 3 presents the errors of the predicted results. For average point-to-curve distance error, the mean value on test set is 1.52 mm. The maximum value of distal tip distance error is 12.33 mm which is only 2.16% of the overall insertion length L . The maximum value of the accumulated angle error is 2.09 degree.

Figure 8 showcases some of the testing results. One-bend shape #54 had the maximum average point-to-curve distance error and #70 had the minimum average point-to-curve distance error. The predicted shapes in Figure 8A and the predicted curvatures in Figure 8B matched the actual shape well.

The 5-fold cross-validation results are presented in Table 4. Again, the training and testing scores are close to each other, indicating that the proposed network worked equally well on all data sets. Note that Dataset 1 contains only shapes with positive curvatures, and the network might have only learned the magnitude correlation between the insertion force and curvature.

We also trained another neural network which has a two-layer unidirectional LSTM (forward only) with MLP of two hidden layers. Selected results Shape #24 and shape #99 are illustrated in Figure 9 for an example. It can be observed that the prediction before the peak is biased and fluctuated for shape #24 and #99 (the blue and red arrows). Furthermore, as the insertion length increases, the longer the length LSTM could use, the better the prediction is. However, the error accumulated in the beginning part of curvature prediction diverges the reconstruction based on it. The mean and max point-to-curve distance is 8.92 and 46.96 mm. The mean and max accumulated angle error is 3.09 and 7.15 degree. These errors are much larger than the proposed model.

Based on above results, we can summarize that the proposed model is good at learning the magnitude relationship between the insertion force and curvature.

The errors for two-bend shapes were huge, comparing with one-bend shape. The mean values for the three errors were: 9.03, 36.76 mm and 31.13 degree. We picked #147 and #156 as examples and illustrated in Figure 10. These two shapes had same first bend and opposite second bend. Shape #147 had two positive curvature bends and shape #156 had a positive bend followed by a negative bend.

Then, the predictions of shape #156 were accurate. The predictions of #147 were accurate only at the first bend section. However, in Dataset 2, 32 two-bend data were mixed with 120 one-bend data, which means the training data set was biased with positive curvature sample. The prediction should be accurate on #147 which had two positive curvature bends.

Thus, the correlation of magnitude learned from one-bend shapes could be learned and the curvature of the first bend could be accurately predicted. It also shows that the signs of curvatures cannot be predicted in this wok.

First, a ninth-order polynomial method is deployed for the curve fitting, which causes fluctuations at the channel’s straight part (a straight line). The fluctuation affects the curvature of tube shape calculation, which are labels feed to neural nets. At the same time, the insertion forces as training inputs are smooth and small over this range. This introduced defects in training process.

Second, we had unevenly distributed datasets. As introduced in Section Experimental Setup, the glass fiber reinforced plastic tube was deformed and fixed on the perforated board by cable ties to create the curved channel. The shapes formed depended on the location of the cable ties. The holes on the perforated board result in limited number of samples, and leading to uneven distribution of the curvature along the channel κ ( s ) . Moreover, for two bend shapes, we had few choices on the angle of first bend section and the curvature of second bend had different signs. These also led to an unevenly distributed dataset.

Third, the size of the data sets was small and the shapes were simple. We measured 120 one-bend shapes, with angles ranging from 0 to 87 degree and curvature ranging from 0 to 6.8 m⁻¹. The range and samples were limited. The one-bend and two-bend shapes were really simple. Moreover, when we attempted the two-bend shapes, the shape of the second bend cannot be accurately predicted for all cases, which might be solved by increasing the amount and range of the data.