Finite-Time Contraction Control of a Ring-Shaped Soft Pneumatic Actuator Mimicking Gastric Pathologic Motility Conditions

Abstract

Soft gastric simulators are the latest gastric models designed to imitate gastrointestinal (GI) functions in actual physiological conditions. They are used in in vitro tests for examining the drug and food behaviors in the GI tract. As the main motility function of the GI tract, the peristalsis can be altered in some gastric disorders, for example, by being delayed or accelerated. To simulate the stomach motility, a GI simulator must achieve a prescribed healthy or pathological peristalsis. This requires the simulator to be controlled in a closed loop. Unlike conventional controllers that stabilize a controlled plant asymptotically, a finite-time controller regulates state variables to their equilibrium points in a predetermined time interval. This article presents the design and implementation of a finite-time, model-based state feedback controller (based on the differential Riccati equation) on a soft robotic gastric simulator's actuators for the first time. We propose a mass-spring-damper model of a ring-shaped soft pneumatic actuator (RiSPA). RiSPA is a bellows-driven, elastomer-based actuator developed to reproduce motility functions of the lower part of the stomach (pyloric antrum). The proposed model is augmented by a new approach for modeling the soft tissues, where the moments of inertia of the system constituents are considered as time-varying functions. The finite-time controller is successfully applied on the RiSPA in numerical simulation and experimental implementation, and the results were thoroughly analyzed and discussed. Its accuracy and the ability to control in a predetermined time are highlighted in the tracking of peristalsis trajectory and contractive regulations.

Introduction

Food and pharmaceutical industries are demanding gastrointestinal (GI) simulators, because they are cheaper and more efficient on a large scale compared with in vivo tests.^1,2 Soft robotic modules and actuators are used to develop the new generation of GI simulators with a near-realistic simulation ability to mimic GI tract functions.^1–5 However, due to some difficulties like lacking the rigid skeleton, complex contractive movements, and embodied sensory system, soft robotic GI simulators became challenging in the design, fabrication, and control.

The primary motility function in the GI tract (particularly in the stomach) is called the peristalsis, an advancing wavelike motion that breaks down food particles by continuous radial contractions and relaxations of the muscles. To reproduce a normal peristaltic wave, various actuation strategies like rope driven, syringe driven, or air chamber were introduced for soft robotic stomach simulators, all of which used open-loop control approaches.^1–3 However, the stomach's peristalsis and contractive movements have specific dimensions,^6,7 which require a feedback system to be reproduced precisely by stomach simulators.

In addition, certain diseases lead to indigestion and disturb the digestion motility process, urging stomach muscles to work with abnormality. Rapid gastric emptying (dumping syndrome) is a stomach abnormality that usually happens to people who have certain types of gastric surgery, leading to the loss of duodenal feedback.⁸ Dumping syndrome accelerates the stomach motility and makes the stomach empty its content into the small intestine faster than normal. Gastroparesis (partial paralysis of the stomach) is another disease common in people with diabetes. By damaging the vagus nerve, which is responsible for the contractive movement, gastroparesis prevents stomach muscles from working normally.⁹

Another gastric disorder, called functional dyspepsia, disturbs peristalsis during and after a meal. These abnormalities delay or accelerate gastric emptying, induce abnormal antral contractions, and accommodation (the ability of the stomach to spread out appropriately to the size and timing of a meal) issues in the fundus and pyloric antrum.¹⁰ Besides, during the digestion process, the stomach is exposed to a range of disturbances from the food contents, which also demands implementing a closed-loop control system for a stomach simulator. Therefore, if both the healthy and pathologic stomach motilities are of interest, the stomach simulator must be able to generate contractive movements in any fixed time interval.

This can be interpreted as the regulation or tracking problem in finite time in the control perspective. While conventional closed-loop controllers do not guarantee that the system gets regulated to its final (desired) states during a specific time, finite-time controllers can regulate the system in any fixed time interval. Should contractive movements be accelerated or delayed (compared to the normal circumstance), the controller must generate input signals acting in shorter or longer time intervals, respectively.

Finite-time control schemes can stabilize the dynamic system in any finite-time interval, while the state variables remain in a limited bound. Amato et al. solved finite-time output and state feedback problems for systems with uncertainties and unknown disturbances and proposed sufficient conditions for finite-time stability using linear matrix inequalities (LMIs).^11–13 Combining the H-infinity method and LMIs, Liu et al. proposed a finite-time control and applied it on a hydraulic turbine governing system model.¹⁴ They showed that the finite-time control has a concise transient state and a stable time. Using the differential Riccati equation, the finite-time H-infinity output controller was designed for continuous linear dynamic systems.¹⁵

The state dependent Riccati equation (SDRE), as a promising model-based feedback control strategy, allows the nonlinearities in the model and possesses the ability to act as a finite-time controller.¹⁶ The conventional SDRE turns the nonlinear structure of the model into a linear form and minimizes a quadratic, infinite horizon performance index to achieve the control law through solving the algebraic Riccati equation. To have a performance in finite time, the conventional infinite horizon performance index is replaced by a finite horizon index along with a penalty term (representing state variable vector in the final time), and the algebraic Riccati equation is converted into a differential equation.^17–21 To design and implement the SDRE control, developing the dynamic model is the first step.

Dynamic models for many soft manipulators and actuators (including ring-shaped actuator models to mimic contractive movements) were developed using first principles.^22–27 This work develops a nonlinear dynamic model of a bellows-driven, elastomer-based ring-shaped soft pneumatic actuator (RiSPA). RiSPA was designed and fabricated as the primary actuation system for stomach simulators, reproducing the contractive movements of the pyloric antrum.²⁷ The proposed dynamic model captures all the actuator's behaviors, along with the compliant nature of the elastomer, and can be readily converted into the state-space form. This model makes it possible to apply a wide range of model-based control schemes on the RiSPA. Then, to replicate the contractive motion of stomach muscles in both the healthy and pathologic modes, the finite horizon SDRE (FH-SDRE) strategy is designed and implemented on the RiSPA. This finite-time controller enables the RiSPA to imitate the peristalsis in a wide range of paces with high accuracy.

The remainder of this article is organized as follows. The nonlinear, mass-spring-damper (MSD) based model of the RiSPA is proposed in The Dynamic Model Development section. In this model, in addition to other deformations, the variations of the moments of inertia of the system constituents are considered. In addition, the model identification is done based on the gray box approach. In FH-SDRE Control section, the FH-SDRE control is designed for the RiSPA based on the proposed nonlinear model. The experimental setup and numerical simulation structures of the FH-SDRE are explained in detail and successfully applied on the RiSPA in Results and Discussions section. This section is primarily divided into tracking the peristalsis trajectory and the regulation problems. The experimental and simulation results are thoroughly discussed and analyzed. Conclusion and Future Works section concludes the article with future works.

The Dynamic Model Development

The anatomy of the stomach is categorized into four main parts, including the fundus and the body (to store the produced digestion gas and food content), the pyloric antrum (the main part regarding the motility function), and the pyloric sphincter (the lowest part of the stomach) (Fig. 1). While the pyloric antrum generates contractions and propels the food content back and forth, the pyloric sphincter restricts the volume and size of the digesta that can enter the intestine. Figure 1 illustrates how the pyloric antrum was estimated by a physical soft actuator called RiSPA.²⁷ The RiSPA comprises five silicon bellows (air springs), as the primary actuation system, enclosed in a polydimethylsiloxane frame. When the pneumatic pressure applies on the bellows, the air pressurizes its hollow chamber, leading to symmetrical expansion, mostly in the longitudinal axes. A soft layer connects all five bellows (made of Ecoflex 0030), and a range sensor is installed at the bottom of each bellows to measure the top part's displacement (Fig. 1).

FIG. 1.

The stomach schematic and its four different main parts. The red-dashed line is the imaginary center line of the stomach. The pyloric antrum is estimated by a ring-shaped actuator with five embedded bellows, with the thickness, inner diameter, and outer diameter of 20, 52, and 82 mm, respectively. A cylindrical soft layer (the green cylinder layer) connects all five bellows, with a diameter of 50 mm, a width of 30 mm, and a thickness of 1 mm.

In previous studies,^25,26 the stomach contractions were not assumed lively enough to be considered a dynamic system and categorized as a quasi-static model (with negligible dynamic effects). In quasi-static models, the time dependency is ignored, for such models run really deliberate. However, in the actual stomach, contractions are rapid and time dependent.^2,28 In addition, some diseases make the stomach movements slower or accelerated.^8–10 Thus, omitting inertia effects and time dependency in the dynamic model makes it unrealistic and incomplete. The finite-time ability, increasing the accuracy, significant reduction in retuning cost and time, the ability to track the undetermined real-time trajectories, dealing with unpredictable changes in the system, and generally achieving better performance demand for designing a finite-time, model-based control, which requires deriving a dynamic model for the system.

Based on the fabricated RiSPA (Fig. 2A), in this article a dynamic model was developed, in which time-dependent movements and deformations were captured. The RiSPA was broken into five 2D hollow circular sectors (HCSs), which were actuated individually. To model the bellows’ deformations, a set of translational springs and dampers was added to the top of each HCS (Fig. 2B). In addition to the radial movements of each HCS ( $r_{i} (t)$ ), the lateral movements were considered ( $θ_{i} (t)$ ). To address the lateral movements, a set of angular (torsional) springs and dampers was added to the bottom of each HCS as well. By importing the lateral movements (), the moments of the inertia of each HCS appeared in the dynamic model subsequently.

FIG. 2.

(A) The prototype of the RiSPA and the schematic of an HCS. (B) The entire model of the RiSPA, including the schematic model of each bellows, as the primary actuation system. HCS, hollow circular sector; RiSPA, ring-shaped soft pneumatic actuator.

Although translational and torsional sets of spring/damper imitate the movements of the RiSPA, they do not accurately represent the time-variant geometry and deformability of the compliant material. Addressing this issue, the HCSs’ moments of inertia, as a function of time, were considered in the model, leading to high nonlinearity of the model. Note that the model is considered decoupled, because the soft layer between each two adjacent HCSs is considerably softer than the bellows and has only 1 mm thickness.

To derive the dynamic model of the RiSPA, the Euler–Lagrange method was used as depicted in Equation (1), in which the generalized coordinate vector was defined as are radial and lateral movements of the ith-HCS, respectively (Fig. 2B). The pneumatic forces applied on each bellows chamber can be presented as $u_{i} (t), (i = 1, \dots, 5)$ , which is the work done on the ith-HCS by an external input force in the direction of $r_{i} (t), (i = 1, \dots, 5)$ (the blue arrows in Fig. 2B). The RiSPA's Euler–Lagrange method is given by the equations $\begin{matrix} \frac{d}{d t} (\frac{\partial L}{\partial {r_{1}}^{'} (t)}) - \frac{\partial L}{\partial r_{1} (t)} = u_{1} (t), \frac{d}{d t} (\frac{\partial L}{\partial {θ'}_{1} (t)}) - \frac{\partial L}{\partial θ_{1} (t)} = 0, \\ \begin{matrix} . \\ . \\ . \end{matrix} \begin{matrix} . \\ . \\ . \end{matrix} \\ \frac{d}{d t} (\frac{\partial L}{\partial {r'}_{5} (t)}) - \frac{\partial L}{\partial r_{5} (t)} = u_{5} (t), \frac{d}{d t} (\frac{\partial L}{\partial {θ'}_{5} (t)}) - \frac{\partial L}{\partial θ_{5} (t)} = 0, \end{matrix}$ (1)

where L is the Lagrangian function, which is the subtraction of the potential energy (U_i) from the kinetic energy (T_i): $\begin{matrix} L_{r} & = \sum_{i = 1}^{5} T_{i} (r_{i} (t), {r'}_{i} (t)) - \sum_{i = 1}^{5} U_{i} (r_{i} (t)); (i = 1, \dots, 5), \\ L_{θ} & = \sum_{i = 1}^{5} T_{i} (θ_{i} (t), {θ'}_{i} (t)) - \sum_{i = 1}^{5} U_{i} (θ_{i} (t)); (i = 1, \dots, 5), \\ L & = L_{r} + L_{θ} . \end{matrix}$ (2)

Consequently, by combining Equations (1) and (2), and choosing the state variable vector as $x (t) = {[\begin{matrix} q (t) & q' (t) \end{matrix}]}^{T}$ and the input vector as , the input-affine, nonlinear dynamic model (of the order of 20) was obtained in the state-space form of $\begin{matrix} x' (t) = f (x (t)) + g (x (t), u (t)) = \\ [\begin{matrix} q' (t) \\ - M^{- 1} (x (t)) [K q (t) + C q' (t)] \end{matrix}] + [\begin{matrix} 0_{10 \times 1} \\ - M^{- 1} (x (t)) u (t) \end{matrix}], \end{matrix}$ (3)

in which: $\begin{matrix} M (x (t)) & = [\begin{matrix} m & 0_{5 \times 5} \\ 0_{5 \times 5} & I (x (t)) \end{matrix}], K = [\begin{matrix} k_{t} & 0_{5 \times 5} \\ 0_{5 \times 5} & k_{θ} \end{matrix}], \\ C & = [\begin{matrix} c_{t} & 0_{5 \times 5} \\ 0_{5 \times 5} & c_{θ} \end{matrix}] . \end{matrix}$ (4)

In Equation (4), $M (x (t)), K, a n d C$ denote the inertia, stiffness, and the damping matrices of the dynamic model, where $m, k_{t}, k_{θ}, c_{t}, c_{θ}$ are diagonal time-invariant matrices of mass, stiffness, and damping coefficient, respectively. On the contrary, $I (x (t)) = d i a g (I_{p} (r_{1} (t)), I_{p} (r_{2} (t)), I_{p} (r_{3} (t)), I_{p} (r_{4} (t)), I_{p} (r_{5} (t)))$ is a time-dependent matrix, including nonlinear, time-variant scalars of , representing inertia and geometry variations of HCSs: $\begin{matrix} I_{x} (r_{i} (t)) = \frac{m g}{4} (R^{2} + {r_{i}}^{2} (t)) (1 - \frac{sin α cos α}{α}), \\ I_{y} (r_{i} (t)) = \frac{m g}{4} (R^{2} + {r_{i}}^{2} (t)) (1 + \frac{sin α cos α}{α}) \\ - \frac{2 m g (R^{3} - {r_{i}}^{3} (t)) sin α}{3 α (R^{2} - {r_{i}}^{2} (t))}, \\ I_{p} (r_{i} (t)) = I_{x} (r_{i} (t)) + I_{y} (r_{i} (t)) \\ = m [0.5 (R^{2} + {r_{i}}^{2} (t)) - 0.38 (\frac{R^{3} - {r_{i}}^{3} (t)}{R^{2} - {r_{i}}^{2} (t)})^{2}] . \end{matrix}$ (5)

The RiSPA's specifications are given in Table 1, where its values were obtained from either the physical measurements (dimensions of the bellows and frame) or the model identification (the mass, the stiffnesses, and the damping coefficients). To characterize the model, a unit step function was fed to the RiSPA, and the step response was captured for each HCS individually (Fig. 3). The step response (over 5 s) for each HCS was repeated several times, while the standard deviation turned out negligible.²⁹

FIG. 3.

Step response analysis of the RiSPA. The figure includes both the experimental and analytical data for specification of the model. The first five legends represent the data from the experimental step response of the RiSPA. The last two legends are from the analytical solution of the dynamic model, which were tuned until approximately matched the experimental data.

Table 1.

The Physical Characteristic for the Ring-Shaped Soft Pneumatic Actuator

Description	Parameter	Value	Unit
Translational stiffness	k_t	0.065	N·mm⁻¹
Angular stiffness	$k_{θ}$	0.068	N·mm·rad⁻¹
Translational damping coefficient	c_t	1	N·s·mm⁻¹
Angular damping coefficient	$c_{θ}$	1.9	N·s·rad⁻¹
Mass of ith-HCS	m	1	gr
Length of bellows	B_l	30	mm
Diameter of bellows	B_d	10	mm
Width of the frame	L	15	mm
Fixed radius of the frame	R	25	mm

HCS, hollow circular sector.

The lack of information in the experimental data can be substituted with the knowledge about the dynamic model. Thus, the analytical response of the model represented in Equation (3) was calculated, while some of the parameters were fixed due to the physical measurements. The combination of analytical information and experimental data (gray box method) can improve the model's accuracy where prior knowledge of the system is not comprehensive enough for satisfactory modeling. The rest of the model parameters were tuned using the MATLAB curve fitting function (cftool) until the analytical response approximately matched the experimental response (Fig. 3). Consequently, by inserting the values from the gray box method and the physical measurements into the dynamic model [Eq. (3)], the nonlinear, time-variant model of the RiSPA, resembling the pyloric antrum, was obtained.

FH-SDRE Control

The FH-SDRE controller, as a model-based control method, was designed based on the proposed nonlinear model of the RiSPA in Equation (3) to stabilize its state variables (radial and lateral movements of each HCS) concerning the regulation and tracking trajectory problems. Since the RiSPA is a continuous physical system, represented by a set of ordinary differential equations [Eq. (1)], its vector functions, including $f (x (t))$ and $g (x (t), u (t))$ , are piecewise continuous and satisfy the Lipschitz condition. Therefore, according to the optimal control theory,³⁰ the control law (input vector) had to minimize the cost function in $t \in [\begin{matrix} 0 & t_{f} \end{matrix}]$ $J (x (t), u (t)) = \frac{1}{2} (\begin{matrix} x^{T} (t_{f}) S x (t_{f}) + \\ \int_{0}^{t_{f}} (x^{T} (t) Q_{c} x (t) + u^{T} (t) R_{c} u (t)) d t \end{matrix}) .$ (6)

In Equation (6) the state cost matrices $Q_{c}$ and $S$ [for regulating the state variables $x (t)$ ] must be symmetric non-negative, and the input cost matrix $R_{c}$ [for regulating the control law $u (t)$ ] must be symmetric positive definite.¹⁶ The quadratic scalar cost function in Equation (6) compromises between the control aim (converging the error of HCSs’ displacements to zero) and the control effort (which is interpreted as the equivalent consumed voltage for pneumatic valves to generate the pneumatic forces during the control task). By factorizing Equation (3) with respect to $x (t) a n d u (t),$ so-called the state dependent coefficient (SDC) parameterization, the nonunique state-space representation of the RiSPA, was obtained as $\begin{matrix} x' (t) = A (x (t)) x (t) + B (x (t)) u (t), \\ A (x (t)) = [\begin{matrix} 0_{10 \times 10} & I_{10 \times 10} \\ - M^{- 1} (x (t)) K & - M^{- 1} (x (t)) C \end{matrix}], \\ B (x (t)) = [\begin{matrix} 0_{10 \times 10} \\ M^{- 1} (x (t)) \end{matrix}], \end{matrix}$ (7)

where the pair of is pointwise controllable for all By constructing the Hamiltonian, including the multiplier (costate variable),

and using the optimally necessary condition,

the control law which minimizes the cost function of Equation (6) was derived as

By substituting Equation (10) into the Hamiltonian and combining the Hamilton-Jacobi-Bellman equation, the following differential equation was obtained: $\begin{matrix} - P' (x (t), t) = & P (x (t), t) A (x (t)) + A^{T} (x (t)) P (x (t), t) \\ - P (x (t), t) B (x (t)) R_{c}^{- 1} B^{T} (x (t)) P (x (t), t) + Q_{c} . \end{matrix}$ (11)

To solve Equation (11), which is known as the differential Riccati equation, and find the control gain $P (x (t), t)$ , a handful of mathematical practices were developed.¹⁶ In this article, the state transition matrix method was used.³⁰ To solve the differential Riccati equation, the following equation was defined based on the Hamiltonian and the costate variable, λ(t)as

which is an ordinary matrix differential equation that can be simply solved by eigenvalues and eigenvectors (leading to an exponential response). Besides, the solution for Equation (12) would have the form of

From the optimal boundary condition (t_f free and x(t_f) fixed),³⁰ and according to the cost function represented in Equation (6)

Given the fact that the costate variable is a function of $P (x (t), t)$ and by substituting Equations (14) into (13), the equation

holds, which gives a numerical solution for the differential Riccati equation. To prove the stability, choosing $V (x (t), t) = x^{T} (t) P (x (t), t) x (t)$ as the positive Lyapunov candidate, the stability of the system was guaranteed in the sense of Lyapunov.

Results and Discussions

Simulation and implementation architecture

To investigate and analyze the efficiency and practicability of the proposed nonlinear dynamic model and the FH-SDRE controller, the numerical simulation (MATLAB programming) was conducted. The cost matrices [ $Q_{c}, S, a n d R_{c}$ from Eq. (6)] were chosen through the numerical simulation safely and efficiently. According to the flowchart in Figure 4A, after defining the required vectors and matrices [Eqs. (3) and (4)], and practicing the SDC parameterization, the cost matrices were added to the system. Since the system was assumed full-state feedback, the main diagonal of $Q_{c}$ was chosen to be nonzero values to cover all state variables. The cost matrices were tuned manually until the results satisfied the control aims.

FIG. 4.

The numerical simulation and experimental setup architecture. (A) The flowchart for the numeric simulation of the FH-SDRE method in the MATLAB. (B) The connection between different hardware and software parts in the experimental setup. For the pressurized state, electro-pneumatic regulators (ITV0030, SMC) and for the vacuum state electronic vacuum regulators (ITV0090, SMC) are used. The relationship between the input signal and the pressure/vacuum outputs is linear in both ITV0030 and ITV0090. The pressure range for the vacuum and pressure regulators is −100 and 500 kPa, respectively. A micro-3-port solenoid valve (005 series, KOGANEI) is used as the switch between the pressure and vacuum states. (C) The block diagram for implementing the experimental tests. FH-SDRE, finite horizon state dependent Riccati equation.

Should the cost of the process (the control effort) outweigh its time duration (the pace of the peristalsis), $R_{c}$ and $Q_{c}$ must be increased and decreased, respectively. In this case, the RiSPA actuates slowly while consuming the minimum control effort. In contrast, if a precise regulation (converging to the desired values in a specific, predetermined time duration) was of importance, regardless of the expenditure of control effort (which means reducing $R_{c}$ ), $Q_{c}$ must increase. The following cost matrices were tuned through the simulation and were used for entire simulations and experimental tests:

The FH-SDRE control was experimentally implemented on the RiSPA to verify the model and the applicability of the finite-time controller. Figure 4B depicts the various sections and connections among data processing (software) and the physical setup (hardware) for the experimental implementation. The control process initiates in the LabVIEW interface. The control signals (see the control unit in Fig. 4C) are sent to the pneumatic valves through the MyRIO board. Then the pneumatic valves convert the control signal (voltage) into air pressure to actuate the RiSPA. By RiSPA's actuation, the sensory system (range sensor (VCNL4040), Fig. 4B) sends a feedback signal to the LabVIEW (PC) through MyRIO and generates the displacement error. Subsequently, the error is fed to the controller and closed the control loop. To protect the air chambers from rupturing under high pressure, we added a saturation bound before feeding the signal into the pneumatic valves.

The peristalsis trajectory tracking

The RiSPA was used to track a peristalsis trajectory that is defined according to physiological features of the pyloric antrum motility. Given HCSs’ initial position of $r_{0} = 6 m m$ (when no pressure or vacuum exists inside the bellows chamber), and the frequency of the antral contraction wave of 0.05 Hz (three cycles per minute),²⁸ the peristalsis trajectory was specified as follows, $\begin{matrix} r_{D} (t) = r_{0} + 4 \times sin (\frac{π}{10} t), \\ {r'}_{D} (t) = 4 \times \frac{π}{10} cos (\frac{π}{10} t), \end{matrix}$ (17)

where $r_{D} (t)$ and $r'_{D} (t)$ are the trajectory and its derivative, respectively. As shown in Figure 5, the control system favorably tracks a peristalsis trajectory with a negligible error. Figure 5 illustrates that the FH-SDRE forced the ith-HCS to follow the peristalsis in both the position and velocity (rate) trajectories. The initial position of the HCS was chosen 12 mm away from the initial position at $r_{0} = 6 m m,$ to show the ability of the FH-SDRE to track the trajectory at the beginning. Figure 5B illustrates that the ith-HCS reached the peristalsis trajectory in <2 s. The root mean square error (RMSE) between HCSs’ displacements (under the FH-SDRE control) and the trajectory, the difference between the blue and the red dashed lines in Figure 5 was calculated as 0.07 mm, which is insignificant against the model's dimensions given in Table 1.

FIG. 5.

The numerical simulation results for tracking a peristalsis trajectory for ith-HCS. (A) The FH-SDRE follows the peristalsis trajectory closely. (C) The FH-SDRE follows the derivative of the trajectory (rate of the displacement [velocity]). (B, D) The magnified illustrations show how the controller tracks and follows the desired trajectory in the beginning.

To verify the tracking trajectory simulation results and the proposed nonlinear model, an experiment was implemented on the RiSPA with the same peristalsis trajectory [Eq. (17)] and initial condition (Fig. 6). The average RMSE among five HCSs for tracking the peristalsis trajectory was calculated as 1.1 mm, which is 10.3% greater than the error from the simulation in Figure 5. Its main reason is the fluctuations that occurred during the experimental tracking trajectory (Fig. 6A). During the implementation of the FH-SDRE, there was a quick rest between each two-consecutive feedback (20 ms), in which no pressure applies.

FIG. 6.

The experiment is implemented following the structure depicted in Figure 4. (A) The experimental results for five HCSs of the RiSPA and how they track the peristalsis trajectory. (B) The average of control efforts (energy consumption) of five HCSs of the RiSPA, required for following the peristalsis trajectory. Note that during the implementation in the LabVIEW programming, the range of the input signal (saturation bound) to protect the actuators was considered 2 V for ITV0030 and −10 V for ITV0090 (the maximum input voltage is ±10 V for both the regulators).

Since the RiSPA's bellows act like a spring, they tend to retreat instantly when the pressure is discharged between those rests, leading to the fluctuations in Figure 6. Another reason for this difference is that the time-variant FH-SDRE gain, $P (x (t), t),$ was generated in real time, based on the feedback from the actual RiSPA, while in the simulation, $P (x (t), t)$ was calculated based on the RiSPA's dynamic model [Eq. (3)]. However, the dynamic model was derived based on first principles [Eqs. (1–5)], and some of its specifications were estimated (Table 1). This creates some differences between the actual RiSPA and the model, leading to the difference in values of $P (x (t), t) .$

The control norm (the magnitude of the input signal vector) was used to measure the control effort. Control norm is a dimensionless criterion representing the cost expenditure during the control process. The average of the control norm during the experimental tests (377.11) was less compared with the simulation (1926.8), which also justifies the difference between their RMSEs.

Although the normal frequency of the antral contraction wave is 0.05 Hz, to validate the performance of the FH-SDRE a comparative experiment on higher velocities was implemented. Having a clear comparison, only the movements of 1st-HCS were captured during different frequencies in 100 s time interval. Figure 7 illustrates FH-SDRE's peristalsis trajectory tracking performance for 0.05 Hz (normal wave), 0.5 Hz (rapid wave), and 1 Hz (very rapid wave). According to Figure 7, the HCS successfully follows the peristalsis under the proposed controller by different frequencies.

FIG. 7.

The comparative experiment test for different frequencies. Three different ranges of frequencies, including 0.05 Hz (20 s period) as the normal speed, 0.5 Hz (2 s period) as the rapid speed, and 1 Hz (1 s period) as the very rapid speed, are compared to each other.

The RMSE was calculated as 1.05, 1.66, and 3.63 mm for 0.05, 0.5, and 1 Hz trajectories, respectively. Due to the actuator's physical limitation and its soft, relatively sluggish nature, the HCS could not retract instantly while tracking the very rapid trajectory, leading to a greater error. Given the average 1.1 mm calculated error and the fact that the entire back and forth stroke of each HCS is 20 mm (Fig. 2A), it can be said that the stomach simulator works averagely with 93% accuracy for 3 to 30 cycles per minute (almost 95% accuracy for normal antral contraction waves). Since the feedback system is of real time, the control system would also deal with variable frequency trajectories within the actuator's physical limitation. This allows us to get real-time trajectories from the actual stomach and feed it to the stomach simulator for mimicking.

Contractive regulation

To show the ability of the FH-SDRE to control the system in finite time, the regulation problem for the contractive movement was solved. Note that cost matrices [Eq. (16)] and physical characteristics (Table 1) were kept the same as the tracking problem. For the simulation, the initial and final state vectors for all HCSs embedded in the RiSPA were chosen $x (0) = {[\begin{matrix} \begin{matrix} R & R & R & R & R \end{matrix} & 0_{1 \times 15} \end{matrix}]}^{T}$ and $x (3) = {[\begin{matrix} 5 & 10 & 17.5 & 15 & 12.5 & 0_{1 \times 15} \end{matrix}]}^{T},$ respectively. The time of the regulation, in which HCSs regulate from $x (0)$ to $x (3)$ , was chosen as 3 s. For the conventional SDRE, which asymptotically stabilizes the system, this time interval is considered short. Figure 8 illustrates how each RiSPA contracts under the FH-SDRE control law.

FIG. 8.

The regulation process for both the conventional SDRE and FH-SDRE in 3 s, where the back-and-forth displacements of HCSs were converted into the radial contraction. (A) The conventional SDRE regulation and (B) The FH-SDRE regulation. To illustrate the controller's ability to regulate to different final states, the main diagonal values for $x (3)$ were chosen differently, leading to an asymmetric contraction. Since the proposed model is decoupled, this does not make a coupling problem among HCSs. According to Table 1, the maximum back-and-forth displacement for each HCS is from R = 25 mm to 0 (100% contraction). Therefore, HCSs initiated to contract from the fully retracted state toward the center of the RiSPA, reaching 80%, 60%, 30%, 40%, and 50% of full contraction (red-dashed lines).

Although in both Figure 8A and B the control laws were the same [Eq. (10)], the control gains were obtained differently in each of them. For the conventional SDRE, the control gain, $P$ , was calculated from the algebraic Riccati equation (time-invariant matrix).¹⁶ However, for the FH-SDRE, $P (x (t), t)$ was calculated from the differential Riccati equation [Eq. (11)] and was a state-dependent, time-varying matrix. While the FH-SDRE met the regulation aim in 3 s, the conventional SDRE could not meet the final state during the same time interval. The average steady-state error of five HCSs from the FH-SDRE was 8.65% of the conventional SDREs. Table 2 shows the comparison between the FH-SDRE and conventional SDRE quantitatively in detail.

Table 2.

The Quantitative Analysis of Each Hollow Circular Sector in the Regulation Problem for Both the Conventional and Finite Horizon State Dependent Riccati Equation Methods in 3 s

ith-HCS	HCS₁	HCS₂	HCS₃	HCS₄	HCS₅	Average of five HCSs
FH-SDRE steady-state error (mm)	0.20	0.13	0.05	0.28	0.44	0.22
Conventional SDRE steady-state error (mm)	2.39	1.66	0.92	3.13	4.59	2.54

FH-SDRE, finite horizon state dependent Riccati equation.

Figure 9A and B illustrates how each state variable, including the angular and translational displacement of ith-HCS in the plane, was regulated to its desired value. A vital criterion in finite-time regulation is the pace of the process. Each state variable must start from its initial position, move toward the final state, and stop there. Otherwise, the system would still have inertia and will not stop moving at the final state. Figure 9C, meanwhile, shows that after finishing the regulation process (in 3 s), the pace of all HCSs turned zero. Figure 9D–F exhibits the reliability test for the regulation problem. To ensure that the FH-SDRE would regulate the RiSPA into any final state, and there is no singularity point while solving the Riccati equation, $x (0) = {[\begin{matrix} R & 0_{1 \times 3} \end{matrix}]}^{T} a n d x (3) = {[\begin{matrix} r_{1} (t) & 0_{1 \times 3} \end{matrix}]}^{T}$ were given to the 1st-HCS as the initial and final states (where $25 m m \leq r_{1} (t) \leq 15 m m$ ). Similar to Figure 9A–C, the state variable was regulated successfully considering different final states in 3 s, while no singularity point occurred.

FIG. 9.

The state variables of the dynamic model of the RiSPA, under the FH-SDRE and the reliability test for the contractive regulation. (A) The angular deformation of HCSs. According to initial and final states, they start from the zero and after the regulation stopped at the zero orientation again. (B) Since final positions of five HCSs were chosen 80%, 60%, 30%, 40%, and 50% of the fully retracted RiSPA, ith-HCS (i = 1 … 5) converged to 5, 10, 12.5, 15, and 17.5 mm, respectively. (C) All the HCSs start moving from zero speed and fully stopped by finishing the regulation process's time interval (in 3 s). (D–F) The reliability test for different final positions for a single HCS.

To highlight the ability to control in finite time, both the FH-SDRE and conventional SDRE methods were experimentally implemented on a single bellows actuator.²⁹ Figure 10A illustrates that by applying the FH-SDRE, the bellows actuator for the first time reached its desired value (20 mm) in just 1 s (the rise time, which is the response's time to reach the desired value, is 1 s). However, by implementing the conventional SDRE, the rise time was around 2 s. It means that with the same initial and final states, and cost matrices, while the steady-state error is almost the same, the finite-time controller acts faster than the conventional controller. Figure 10B shows the comparison between the control inputs of two controllers. The control norm for the FH-SDRE and conventional SDRE was calculated as 669.96 and 519.69, respectively. In fact, the finite-time control spent 1.3 times more energy (control effort) than the conventional SDRE to compensate for a rapid regulation.

FIG. 10.

The experimental regulation problem for a single bellows actuator under FH-SDRE and conventional SDRE controllers (following the implementation algorithm depicted in Fig. 4). Single bellows actuator was developed to mimic the stomach muscles. It was made of the same material as RiSPAs.²⁹ (A) The comparative regulation process for the single bellows actuator. The blue and red dashed lines show the rise time for the FH-SDRE and SDRE, respectively. (B) The input signal generated by FH-SDRE and SDRE controllers.

The experimental results were presented for the contractive regulation of the RiSPA under the FH-SDRE in 3 s (Fig. 11). The initial and final states of HCSs were chosen as $x (0) = {[\begin{matrix} R_{1 \times 5} & 0_{1 \times 15} \end{matrix}]}^{T}$ and $x (3) = {[\begin{matrix} 5_{1 \times 5} & 0_{1 \times 15} \end{matrix}]}^{T},$ respectively. The average steady-state error of HCSs was calculated as 0.08 mm, which is acceptable considering that the radius of pylorus never exceeds 80 mm.³¹ In addition, given that the pylorus only allows the food particles with the size of 1–2 mm diameter or smaller to pass through (to enter the small intestine),^1,2,32 the contractive steady-state error was negligible. Table 3 gives the results of Figure 11, quantitatively. The different values for HCSs in Table 3 (particularly the 5th-HCS) were due to the RiSPA's defects and uncertainties introduced during the fabrication process. However, by looking at the steady-state errors, FH-SDRE as a closed-loop method could deal with this imperfection successfully.

FIG. 11.

The experimental data from the contractive regulation, along with the actual RiSPA, (A) Before and (B) after applying the FH-SDRE control in 3 s. Unlike the tracking trajectory implementation, there is no need to define any delay in the contractive regulation problem. This made the regulation's fluctuation very small against the tracking problem.

Table 3.

The Quantitative Analysis of the Experimental Contractive Regulation Problem of the Ring-Shaped Soft Pneumatic Actuator Under the Finite Horizon State Dependent Riccati Equation

ith-HCS	HCS₁	HCS₂	HCS₃	HCS₄	HCS₅
Control effort (control norm)	228.42	221.74	237.35	251.79	666.08
Steady-state error (mm)	0.0715	0.0428	0.0301	0.0681	0.1952
Rise time (s)	0.17	0.3	0.62	0.47	1.02

Conclusion and Future Works

In this study, we reproduced the peristalsis and contractive movements of the pyloric antrum, the lower part of the stomach, for both the healthy and pathologic states using a RiSPA. The RiSPA was developed with the ability of the contractive movements to be used as the primary actuation mechanism for soft robotic stomach simulators. The FH-SDRE method as a finite-time control scheme was successfully designed and implemented on the RiSPA, providing the entire control system with stability in any predetermined time interval. This makes the RiSPA contract at any desired pace in response to the stomach motility disorders like gastroparesis or functional dyspepsia, potentially in real time or offline. A novel dynamic model was proposed to capture all the deformations and geometry variations of the RiSPA and used in designing the FH-SDRE. In addition to MSD equations, moments of inertia of the RiSPA as time-variant parameters were added to the model to have a near-realistic dynamic behavior.

Experimental tests and numerical simulations were designed and conducted in both the tracking trajectory (where some predetermined state trajectory is to be followed) and regulation (where some steady state is to be maintained) problems to prove the effectiveness and capability of the entire control system in different scenarios. It was shown that the RiSPA follows the peristalsis trajectory successfully under the FH-SDRE control, and its RMSE was insignificant compared to RiSPA's dimensions.

The 1.1 mm RMSE was primarily due to the technical issue regarding the delays (quick rests) while receiving the feedback in implementation structure. The advantage of the FH-SDRE (as a finite-time method) over the conventional SDRE was highlighted in the contractive regulation problem through simulation and experimental results. It was shown that the FH-SDRE could regulate the RiSPA in any predetermined short time interval. In addition, the experimental contractive regulation results were compatible and acceptable compared to the actual pylorus measures and functions. In the future, the FH-SDRE control can be extended and implemented in the whole soft robotic stomach simulator. In addition, the physiological motility data of the pathologic stomach can be fed to the FH-SDRE control loop in real time, and the soft robotic stomach simulator can mimic the motility functions simultaneously. In addition, since the RiSPA is a ring circular soft actuator, it can be used to mimic the heart arrhythmia and irregular intestine and esophagus movements.

Footnotes

Author Disclosure Statement

No competing financial interests exist.

Funding Information

Research was supported by the Riddet Institute, a Centre of Research Excellence, New Zealand.

References

Chen

, Xu

, Fan

, et al. Gastric emptying and morphology of a ‘near real'in vitro human stomach model (RD-IV-HSM). J Food Eng, 2016; 183:1–8.

, Fortner

, Kong

. Development of a Gastric Simulation Model (GSM) incorporating gastric geometry and peristalsis for food digestion study. Food Res Int, 2019; 125:1–13.

Dang

, Liu

, Hashem

, et al. SoGut: a soft robotic gastric simulator. Soft Robot, 2021; 8:273–283.

Bhattacharya

, Ali

, Cheng

, et al. RoSE: a robotic soft esophagus for endoprosthetic stent testing. Soft Robot, 2021; 8:397–415.

Dirven

, Chen

, Xu

, et al. Design and characterization of a peristaltic actuator inspired by esophageal Swallowing. IEEE ASME Trans Mechatron, 2013; 19:1234–1242.

Kwiatek

, Steingoetter

, Pal

, et al. Quantification of distal antral contractile motility in healthy human stomach with magnetic resonance imaging. J Magn Reson Imaging, 2006; 24:1101–1109.

Pal

, Brasseur

, Abrahamsson

. A stomach road or “Magenstrasse” for gastric emptying. J Biomech, 2007; 40:1202–1210.

Ukleja

. Dumping syndrome: pathophysiology and treatment. Nutr Clin Pract, 2005; 20:517–525.

Abell

, Bernstein

, Cutts

, et al. Treatment of gastroparesis: a multidisciplinary clinical review: The American Motility Society Task Force on Gastroparesis (members in alphabetical order). Neurogastroenterol Motil, 2006; 18:263–283.

10.

Brun

, Kuo

. Functional dyspepsia. Therap Adv Gastroenterol, 2010; 3:145–164.

11.

Amato

, Ariola

, Dorato

. Finite-time control of linear systems subject to parametric uncertainties and disturbances. Automatica (Oxf), 2001; 37:1459–1463.

12.

Amato

, Ariola

, Cosentino

. Finite-time stabilization via dynamic output feedback. Automatica (Oxf), 2006; 42:337–342.

13.

Amato

, Ariola

, Cosentino

. Finite-time control of discrete-time linear systems: analysis and design conditions. Automatica (Oxf), 2010; 46:919–924.

14.

Liu

, Wang

, et al. Finite-time H-infinity control of a fractional-order hydraulic turbine governing system. IEEE Access, 2018; 4:57507–57517.

15.

Panteleev

, Yakovleva

. Sufficient conditions for H-infinity control on the finite time interval. J Phys Conf Ser, 20211925:1–8.

16.

Nekoo

. Tutorial and review on the state-dependent Riccati equation. J Appl Nonlinear Dyn, 2019; 8:109–166.

17.

Cheng

, Lewis

, Abu-Khalaf

. Fixed-final-time-constrained optimal control of nonlinear systems using neural network HJB approach. IEEE Trans Neural Netw, 2007; 18:1725–1737.

18.

Çimen

, Banks

. Nonlinear optimal tracking control with application to super-tankers for autopilot design. Automatica (Oxf), 2004; 40:1845–1863.

19.

Jung

, Park

, Kim

, et al. Hardware-in-the-loop simulations of spacecraft attitude synchronization using the state-dependent Riccati equation technique. Adv Space Res, 2013; 51:434–449.

20.

Lin

, Vandewalle

, Liang

. Analytical representation of the state-dependent coefficients in the SDRE/SDDRE scheme for multivariable systems. Automatica (Oxf), 2015; 59:106–111.

21.

Heydari

, Balakrishnan

. Path planning using a novel finite horizon suboptimal controller. J Guid Control Dyn, 2013; 36:1210–1214.

22.

Chen

, Xiong

, Liu

, et al. Fabrication and dynamic modeling of bidirectional bending soft actuator integrated with optical waveguide curvature sensor. Soft Robot, 2019; 6:495–506.

23.

Mustaza

, Elsayed

, Lekakou

, et al. Dynamic modeling of fiber-reinforced soft manipulator: a visco-hyperelastic material-based continuum mechanics approach. Soft Robot, 2019; 6:305–317.

24.

Uppalapati

, Krishnan

. Towards pneumatic spiral grippers: modeling and design considerations. Soft Robot, 2018; 5:695–709.

25.

Dang

, Stommel

, Cheng

, et al. A soft ring-shaped actuator for radial contracting deformation: design and modeling. Soft Robot, 2019; 6:444–454.

26.

Dang

, Devaraj

, Stommel

, et al. Experimental investigation into the dynamics of a radially contracting actuator with embedded sensing capability. Soft Robot, 2020; 7:478–490.

27.

Hashem

, Kazemi

, Stommel

, et al. A biologically inspired, ring-shaped, soft pneumatic actuator for large deformations. Soft Robot, 2021. [Epub ahead of print]; DOI: 10.1089/soro.2021.0013.

28.

Schulze

. Imaging and modelling of digestion in the stomach and the duodenum. Neurogastroenterol Motil, 2006; 18:172–183.

29.

Hashem

, Stommel

, Cheng

, et al. Design and characterization of a bellows-driven soft pneumatic actuator. IEEE ASME Trans Mechatron, 2020; 1–12. [Epub ahead of print]; DOI 10.1109/TMECH.2020.3037643.

30.

Kirk

. Optimal Control Theory: An Introduction. Englewood Cliffs, N.J.: Courier Corporation (Prentice-Hall), 2004.

31.

Salessiotis

. Measurement of the diameter of the pylorus in man: Part I. Experimental project for clinical application. Am J Surg, 1972; 124:331–333.

32.

Kong

, Singh

. A human gastric simulator (HGS) to study food digestion in human stomach. J Food Sci, 2010; 75:627–635.