S OFT robots are constructed from soft and compliant mate- rials and structures. This unique nature endows soft robots with competitive advantages such as intrinsic adaptability to diverse environments, enhanced safety and dexterity during services, and flexibility for various applications. Soft robots exhibit a theoretically infinite number of degrees of freedom (DOFs) and thus can generate highly nonlinear behaviors, which leads to challenges in accurately modeling and controlling their deformation, compliance, and behaviors [1]. As a result, soft robotics theoretical models are normally much more complicated than those for traditional rigid robots. To tackle this major problem, remarkable studies have been performed, ranging from physics/geometry-based models and their approximation and discretization [2], numerical methods such as finite element methods (FEMs) [3], sensing enabled morphological and reservoir computing [4], to pure data-driven methods using machine learning techniques [5].
Consequently, several theoretical models have been developed considering the dynamics principles such as continuum mechanics models, geometrical models, discrete material models, and surrogate models [1]. Continuum mechanics models consider soft bodies as a set of continuum particles; thus, the main idea is to equivalently simplify the representative numbers of DOFs for control inputs for driving the motion [6]. Geometrical models are usually constructed from the soft robots' geometrical shapes after applying specific loads [7], hence can describe the deformed body as a space curve defined by mathematical functions or divide the soft body into a set of circular arc shapes. Discrete element models define the soft continuum body as discrete material elements such as springs, dampers, and masses [8]. Physics-based models can provide accurate descriptions of soft robots' nonlinear behaviors; however, they may struggle to integrate seamlessly with real-time sensing and control. Approximation and discretization can reduce the DOFs of such models but may need assumption and simplification and thus may sacrifice accuracy.
Machine learning/deep learning (ML/DL) techniques can enrich soft robotics with models that are able to deal with highly nonlinear problems while maintaining the necessary accuracy [9]. Such data-driven models are suitable for real-time control; however, they rely on extensive and high-quality training data, which can be costly to collect. In addition, applying ML to soft robot models while considering the underlying physics is still one of the major challenges [10]. This is because extracting all the features from only the given data without considering the underlying physical principles may yield unrealistic predictions, especially for those out of the training range.
Physics-guided deep learning (PGDL) offers a compelling approach to overcome the limitations of using either ML or physical models in complex engineering problems. Unlike traditional ML models, PGDL embeds governing physical laws directly into the learning process, ensuring that the predictions align with fundamental physics. Generally speaking, PGDL can address four main problems, solving differential equations, learning dynamics residuals, discovering governing equations, and forecasting dynamic systems [11]. The main idea is to solve partial differential equations (PDEs) that describe the physics of the system using data from analytical and/or numerical solutions and thus have the uniqueness to integrate the strengths of ML with the robustness and interpretability of physics-based models. By utilizing numerical results obtained from physical models such as FEM for learning, this approach substantially lowers the computational effort required for network training, as it is consistent with the underlying physical laws that govern the dataset. This method enhances the prediction accuracy and reduces the reliance on large datasets by embedding physical constraints directly into the learning process. However, coupling physics into deep learning is not yet well studied for soft robot modeling and control.
This study aims to develop a new PGDL-based computational framework that integrates first principle physics and deep learning as surrogate models for soft robotics. We leveraged physics-informed neural network (PINN) [12] that can solve PDEs using spatial coordinates to model intricate soft robotics and established the physical model using the Navier-Cauchy equation to construct the framework. We constructed three soft finger models with different geometries and behaviors as three testbeds. These three chosen soft fingers have been experimentally validated in previously published work. The results demonstrated high accuracy and efficiency in all three cases, which showed that our model can be generalized over a wide range of soft pneumatic systems. In addition, two training and testing approaches were employed using 5% of the total data for training and both covered and uncovered scenarios for testing. Using both approaches, our proposed PGDL model accurately predicted the base curvature and the deformation, demonstrating its robustness and generalization capabilities across different loading conditions. This study is among the first to integrate first principle physics as the governing equation into deep learning for soft robot modeling.
The remainder of this letter is organized as follows: Section II presents the proposed computational framework. Section III covers the FEA implementation for data generation, and Section IV provides the experimental setup and the results for two different training models followed by discussion and conclusion in Section V.
The proposed PGDL aims to fundamentally employ PDE solving in the learning process to capture soft robots' complex behaviors. Once trained, these surrogate models can replace computationally expensive numerical simulations, significantly reducing computation time and enabling real-time control. This section presents the proposed framework in detail.
Traditional applications of PGDL use PINNs to solve PDEs that describe physical models in terms of spatial and temporal variables as long as these physical models can be intrinsically represented by the corresponding governing equations, auxiliary equations, boundary conditions, and initial conditions. Modified PINNs have been attempted in modeling soft robotics, and some researchers have adapted the PINN framework to address the unique challenges of soft robotics. Notably, Sun et al. [13] proposed the Physics-Informed Recurrent Neural Networks (PIRNN) that used a first-order linear system with a rate-independent hysteresis as the physical model and integrated it with a Recurrent Neural Network (RNN) and used the explicit Euler method to predict the output of this physical model. Then this prediction was used as additional input and output to guide the RNN learning. Liu et al. [14] employed Lagrangian and Hamiltonian formulations for modeling soft robot dynamics into NN, which allows for non-collocated control using the generalized coordinates, momenta, and the control input as the input layer.
To generalize the applications of coupling physics and ML/DL for soft robotics, a clear description of the underlying physics with the temporal and spatial coordinates that align the control (or prediction) points with the physical model's behavior is still desired to be well integrated with the learning process, which is not yet well studied.
Pneumatic soft fingers constructed by pneumatic networks (PneuNets) are flexible and compliant actuators made from elastic or hyperelastic materials with integrated inflatable chambers and channels in certain sequences. These soft fingers are often designed and controlled to produce desired motions that can mimic human fingers for grasping and manipulation when well pressurized. Other designs such as bending joints and bellow structures have also been attempted.
This study used pneumatic soft fingers as testbeds to implement and validate the proposed PGDL framework. Three soft finger models with different geometries, materials, and behaviors were constructed as three testbeds, which were tested and validated by other studies [15], [16], [17]. FEA simulations were also performed in these three studies. The first is a silicon-based PneuNets made by [15]; the second is a 3D-printed soft finger with local bending joints [16]; and the third is a bellow-type soft actuator [17]. Fig. 1 shows the detailed design of the first case with the structure of its PneuNets as an example, and the dimensions of this design are detailed in Table I.
TABLE I PARAMETRIC GEOMETRY OF THE PNEUNET
Soft robots' behaviors can be complex due to their highly flexible and deformable nature. The physics underlying such mechanical behaviors of these soft structures made from elastic or hyperelastic materials can be described using the Navier-Cauchy equation [18]. For fully compliant soft robots such as PneuNets, considering an arbitrary volume, the differential form of the equilibrium equation can be given by:
where ∇• is the divergence operator; σ is the Cauchy stress tensor; F is the body force per unit volume caused by gravity. The Navier-Cauchy equation in equilibrium can illustrate that the internal forces within the body and the external forces must balance out to maintain the state of equilibrium. The divergence of the stress tensor, ∇ • σ, represents the internal forces within a material that result from stress, including the nine components of the normal and shear stresses; the stress tensor, σ, describes how internal forces are distributed within the body, and the body force vector, F, represented the external forces.
The Navier-Cauchy equation can be expressed in terms of the three-dimensional components, the indices i, j take the values 1, 2, 3 for three-dimensional x, y, z, using Einstein summation convention as:
where σ ij,j and F i represent the index notation of the stress tensor and the body force density, respectively; and x j represents the three spatial coordinates.
When relating stress to displacement as a dependent variable in the defined domain, Navier-Cauchy equation formulation is achieved using the stress-strain-displacement relationship. The constitutive equation that relates the stress with the strain (Hooke's law) in index notation is:
where λ and µ are Lame's constants; ε ij is the strain tensor in index notation; and δ ij is the Kronecker delta:
Then, the strain-displacement relation can be given as:
where u is the displacement components in the three spatial coordinates and also i, j take the values 1, 2, 3.
Since both the stress-strain and strain-displacement relations are obtained, the stress can be represented in terms of displacement:
Then, the obtained stress can be written as:
The general formulation for describing deformations of elastic bodies in differential form can be given by substituting ( 6) into (2). The Navier-Cauchy equation is then obtained as a function of the displacement vector u and the body force vector F as:
where Lame's constants, λ and µ, are λ = 10.78 MPa and µ = 0.22 MPa. These values were obtained using the Yeoh model, where µ= 2C 10 with C 10 = 0.11 MPa, and λ = 2µν 1-2ν , with an approximation value for Poisson's ratio ν = 0.49 [19].
In this form, the Navier-Cauchy equation contains two main terms. The first term, (λ + µ)∇(∇ • u), represents the volumetric response of the material, modeling the changes in bulk deformation and understanding how the body expands or contracts under an applied force. The second term, µ∇ 2 u, captures the internal shear deformations that explain how the body bends, twists, or experiences other complex deformations, which can be used along with the first term to describe the soft robots' elastic properties and deformation behaviors. By solving this equation within a physics-guided neural network, the model ensures adherence to physical laws, leading to more accurate and reliable predictions of the robot's behavior under various loading conditions.
In the proposed soft finger manipulations as described in Section II-B, the deformations can be considered in two dimensions (2D), x-z plane. Then the governing equation can be given explicitly by:
where u x and u z are the deformations in x and z directions, respectively; and f x and f z are the body forces per unit volume caused by gravity in x and z directions, respectively. This physical model starts from the first principle and provides a general physical model for fully soft and compliant robotics systems. This study focuses on studying the deformation in equilibrium aiming to provide a foundational understanding of how these structures balance internal stresses and external forces. This will offer critical insights into their stability and performance for future more complex transient behaviors.
Coupling the physical model described in Section II-C into the deep learning process can be illustrated in our proposed PGDL framework as in Fig. 2. Through this framework, surrogate models are built via deep learning that incorporates and integrates the above physics by providing two main parts, the ML part through feedforward neural network (FNN) which is designed for static and kinematic modeling, and the physics-guided part through the physical model described using the Navier-Cauchy equation as the loss function, which enforces and constrains the output deformation from the FNN to satisfy and learn from PDEs in (8) that contain the underlying physics.
Regarding the ML part, an FNN was defined by fully connected layers including the input layer, hidden layers, and output layer to provide a complex multi-dimensional group of nested functions, which perform a linear transformation followed by a nonlinear activation function to capture the non-linearity relationship from the training data:
where h is the NN output layer, which contains the deformations u x in x direction and u z in z direction; A is the NN input layer, which includes the pressure p and the soft finger geometry x and z; R i are the activation functions, where i = 1, 2, . . ., M ; and M is the number of transformational layers where weights and activation functions are applied including hidden and output layers; and W i and b i are the edges' weights and biases, respectively.
The physics-guided part includes two residual functions to be minimized as f (x, z) and g(x, z):
Then, it proceeds by approximating u x and u z utilizing a single NN with two outputs. Incorporating this prior assumption along with (10), it leads to a physics-guided deep learning model [f (x, z)| p , g(x, z)| p ], where f and g are functions of the geometry, conditioned on the specific pressure p that caused the deformation of the selected geometry node x and z. The goal of the training process is to train the NN parameters (u x , u z )| p and (f (x, z)| p , g(x, z)| p ), for each corresponding pressure p by minimizing the summation of the MSE of the total loss function L:
where L u x ,u z and L f,g are the loss terms for the deformations and the residual functions; {p i , x i , z i , u i x , u i z } N i=1 represent the overall training (actual) data as vectors with a size of N =20,000, including the pressure p that caused the deformation of the selected soft finger geometry (x and z); exact deformation in x direction (u i x ) and z direction (u i z ); and ûi x and ûi z are the predicted deformations in the x and z direction, respectively.
The leading architecture of the PGDL framework was built using an FNN. In this static equilibrium problem, the PGDL consists of one input layer, eight hidden layers, and one output layer. The input layer contains the geometry of the soft finger in the x and z planes (x and z) and the applied pressure (p). The hidden layers were tested with different sizes and ended by using eight layers with 100 neurons. The output layer is designed to predict the deformation in the x direction (u x ) and the deformation in the z direction (u z ). Rectified Linear Unit (ReLU) was used as an activation function due to its ability to handle the non-linearity behavior in the data, thus effectively capturing the nonlinear behavior of the soft finger with the Navier-Cauchy equation.
The overall amount of data obtained for the finger face from FEM after applying the pressure in 100 sub-steps was 443,900. A total of 20,000 data points, corresponding to a mere 5%, were randomly selected to cover the overall data, forming the size of the inputs x, z, applied pressure p, and the corresponding output ûx and ûz . The training started by flattening all the input and output data into a vector of size (20,000, 1), then normalized in a range of [0, 1] to ensure stable training. After prediction, the network's outputs (û x and ûz ) are reshaped back to their original ranges before computing the derivatives required for the derivatives layer and the physics-based loss functions. This step is essential for accurately applying the neural network's predictions to the PDEs that govern the system and the optimization functions.
Adam optimizer was used to train the network and update the weights with a learning rate of 0.001, which was selected through preliminary experiments to balance convergence speed and stability. A batch size of 512 was implemented throughout the training process, and the training was conducted over 100,000 epochs, with the loss functions guiding the optimization process.
The FEA of the three pneumatic soft fingers described in Section II-B was achieved by numerical simulation using ANSYS, in which the fingers are considered as a uniform continuum body. The silicon-based PneuNets model [15] was chosen to explain the FEA and the initial position was defined perpendicular to the gravity as a horizontal fixation by setting the boundary condition using a fixed support to the left side of the actuator, as shown in Fig. 3. An air pressure was then applied normal to all the internal walls. Contact among surfaces can be expected due to the material's nonlinear deformation and hyperelastic characteristics, so self-contact interactions between the outer walls were defined as frictionless surface-to-surface contact. The material was selected as a second-order hyperelastic Yeoh model with coefficients C 1 = 0.11 MPa and C 2 = 0.02 MPa with a mesh size of 3 mm, to describe the hyperelastic materials.
The simulation results of the curvature shape and the total deformation under the actuated pressure are shown in Fig. 3. Using FEA method, we created a high-resolution dataset for the training process. By applying air pressure and gravitational force at 100 sub-steps, we obtained deformation data under various loading conditions, resulting in 100 snapshots for the deformation in x and z coordinates. The amount of the face nodes after meshing was 4,439 nodes as (x, z) points for 100 pressure snaps (steps) formulated the 443,900 overall training data (4,439×100) for all the inputs and outputs.
Using the generated dataset from FEA, two different training and testing approaches were considered to check the efficiency of the PGDL framework. The first approach was comprehensive as the training data set was created by randomly selecting 20,000 points, a mere 5% of the total available data, to describe the entire domain of the exact solution domain after applying different pressures for each of the three soft pneumatic actuators (0-50 kPa for the first, 0-45 kPa for the second, and 0-120 kPa for the third) and testing with the remaining data. The second approach, which was more challenging, the first actuator was chosen as a case study to be able to further explore the performance of our framework by training the model at 0-44 kPa, using the same amount of data of 5%, and tested up to 50 kPa to be able to judge the ability of the framework to predict uncovered loading conditions and situations. The uncovered loading conditions were chosen to be 45 kPa, 47 kPa, and 50 kPa. Testing for both approaches was achieved by reshaping the overall data to its original 100 pressure snaps as (4,439, 100) and testing the prediction at different snapshots.
A detailed comparison between the FEA results and the PGDL framework is presented under different loading conditions for the three different actuators, seen and uncovered, to highlight the ability of this computational framework to address the limitations of the FEA. The results show that the PDGL framework can provide an accurate prediction and can also work in real-time control applications.
1) First Approach: The main idea of the first approach was to train the model till the maximum applied pressure using a small amount of the available data leveraging the power of the physical model in the learning process. Fig. 4 shows the deformation results in the x direction with different loading conditions for the three different soft fingers the silicon-based PneuNets, the finger with local bending joint, and bellows-type actuator, respectively. The FEA results of the three conditions are presented in Fig. 4 values of FEA results and PGDL predictions were compared. Take Fig. 4(a) and (d) of the silicon-based PneuNets case for example, the maximum and minimum deformations in FEA are 11.234 mm and -63.426 mm, respectively; whereas those for PGDL prediction are 11.26 mm and -63.43 mm. In the bellows-type actuator case, as shown in Fig. 4(c) and (f), the maximum and deformations in FEA are 5.015 mm and -150.38 mm, respectively; whereas those for PGDL prediction are 5.00 mm and -150.63 mm. The MSE for all three representative conditions were between 10 -3 and 10 -4 mm 2 , showing the prediction accuracy. Likewise, Fig. 5 illustrates the deformation results in z direction for the same three soft fingers. Similar to the x direction results, the PGDL predictions are compared with the FEA results, and the maximum and minimum deformations are calculated and compared. As shown in Fig. 5(a) and (d) that are for the silicon-based PneuNets case, the maximum and minimum deformations in FEA are 0.25 mm and -105.24 mm, respectively whereas those for PGDL prediction are 0.23 mm and -105.16 mm. In Fig. 5(b) and (e) of the finger with local bending joint case, the maximum and minimum deformations in FEA are 0.21 mm and -69.94 mm, respectively; whereas those for PGDL prediction are 0.21 mm and -70.10 mm. The MSE for the entire body is still within 10 -3 and 10 -4 mm 2 .
In the experiments, we selected three soft actuators with different geometries and behaviors aiming to estimate the model's capability of predicting the deformation in different conditions without losing its accuracy.
In addition to the maximum and minimum deformation, we also used the curvature of the silicon-based PneuNets baseline to compare the results from FEA and from the PGDL prediction. The comparison under different applied air pressure conditions is shown in Fig. 6 with the air pressure of 10-50 kPa. The solid blue line is the curvature of the finger baseline in FEA whereas the red dash line is the PGDL-predicted finger baseline. The results show that the predicted baseline from the PGDL aligns with the FEA data. As depicted in Fig. 6, the actual and predicted curves overlap precisely for all tested pressures, demonstrating the model's exceptional predictive accuracy.
To quantitatively estimate the prediction accuracy, the MSE between the exact deformation from FEA and the predicted output from our framework was calculated. Fig. 7 lists the MSE for five conditions of the applied pressures of 10 kPa, 20 kPa, 30 kPa, 40 kPa, and 50 kPa, corresponding to the given cases in Fig. 6.
The training was performed with less than 5% well distributed points over the whole body, and the MSE varied between 10 -3 and 10 -4 mm 2 for x and z directions. The beginning and ending TABLE II MSE OF THE UNCOVERED CONDITIONS (UNIT: MM 2 )
cases (i.e., 10 kPa and 50 kPa) had relatively higher MSE but still in the range of 10 -3 mm 2 . The cases in between had lower MSE. The results from the first approach show that using the Navier-Cauchy equation as the physical model to capture the deformation of the soft robotic system and coupling it into the learning process enhanced the model's performance to a large extent for all the loading conditions within the training data range.
2) Second Approach: The second approach was performed using the silicon-based PneuNets soft actuator to perform training with a limited range of loading conditions and predict uncovered air pressure conditions to be able to further explore the full capability, robustness and generalizability of our framework. Specifically, the training was covered up to 44 kPa whereas the PGDL prediction was achieved up to 50 kPa. By training the model with a subset of the data and testing it on higher, uncovered pressures of 45 kPa, 47 kPa, and 50 kPa, we hoped to test the boundaries of the model's predictive capabilities. The results showed that this model was able to predict these uncovered conditions from limited training with low MSE. The experiment and results are detailed below.
The training for the second approach was performed up to 44 kPa and the performance at the three uncovered loading conditions (45 kPa, 47 kPa, and 50 kPa) was tested. The prediction performance was evaluated in terms of the base curvature as shown in Fig. 8. Despite not encountering these specific pressures during training, the model could still achieve high accuracy with the low MSE in the order of 10 -3 mm 2 . Table II lists the MSE values between the exact and predicted deformations for these three uncovered conditions. With such low MSE, the base curvatures for these three pressure conditions were nearly identical to the actual data. Fig. 9 presents the predicted deformation results and the maximum and minimum deformation values at 50 kPa in both x and z directions. The x direction prediction had the MSE of 0.017 mm 2 whereas the z direction prediction had the MSE of 0.007 mm 2 . Compared with the FEA results and PGDL predictions from the first approach (Figs. 4(a) and (d) and 5(a) and (d)), the second approach can achieve comparable predicting performance using the proposed PGDL, which underscores the model's capability of predicting previously unencountered loading scenarios. The consistently low MSE values across these conditions demonstrate the model's robustness and capacity to maintain accuracy even when extrapolating beyond the training data.
This study proposed a physics-guided deep learning framework for soft robot modeling that integrates first principle physics into neural network, which offers a new approach for modeling soft robots' complex, nonlinear behaviors. This study demonstrated the effectiveness of the PGDL framework in accurately predicting the deformation and the base curvature of a pneumatic soft finger across various loading conditions. An FEA model was created to generate a high-resolution dataset by capturing the deformations under different loading scenarios. Two training and testing approaches were utilized the MSE of both approaches were between 10 -3 and 10 -4 mm 2 . This study is among the first to couple first principle physics into deep learning for soft robot modeling, also there is still future work to further improve the PGDL framework.
This study starts with the formulation of equilibrium for predictions and therefore can be used for the prediction of equilibrium conditions at this stage. As a starting step, it provides insights into how soft structures balance internal stresses and external forces and is fundamental to the next-step complex transient behavior modeling using PGDL. This method intrinsically considers physics and has the capability of considering more complex situations by coupling the corresponding governing equations, auxiliary equations, initial, and boundary conditions into deep learning. Therefore, the framework's structure needs to be adaptive to consider these equations and still provide high accuracy. For future work, we will improve this proposed PGDL framework by using the transient model and including the temporal coordinates to capture the behavior over the time domain. This proposed PGDL has a high potential to be applicable and extended to different types of soft robots for soft robot modeling.
Authorized licensed use limited to: University of Virginia Libraries. Downloaded on February 07,2025 at 23:18:52 UTC from IEEE Xplore. Restrictions apply.