A New Hybrid Viscoelastic Soft Tissue Model based on Meshless Method for Haptic Surgical Simulation.

This paper proposes a hybrid soft tissue model that consists of a multilayer structure and many spheres for surgical simulation system based on meshless. To improve accuracy of the model, tension is added to the three-parameter viscoelastic structure that connects the two spheres. By using haptic device, the three-parameter viscoelastic model (TPM) produces accurate deformationand also has better stress-strain, stress relaxation and creep properties. Stress relaxation and creep formulas have been obtained by mathematical formula derivation. Comparing with the experimental results of the real pig liver which were reported by Evren et al. and Amy et al., the curve lines of stress-strain, stress relaxation and creep of TPM are close to the experimental data of the real liver. Simulated results show that TPM has better real-time, stability and accuracy.

INTRODUCTION

Virtual surgery simulation enables physicians to interact with flexibly customized simulate environments based on visual, auditory and haptic feedback without potentially harmful contact with real patients. With the use of virtual surgery, surgeons have the opportunity to surgically train in a low-cost and ethically sound environment [1]. On the other hand, real-time surgical simulation and medical training must employ the physical modeling and behavioral modeling of soft tissue so as to achieve realistic visual and tactile effects [2].

Deformable models can be classified into two categories: physics-based and geometry-based. Physics-based methods are based on continuum mechanics, and can get accurate simulation results by directly solving the partial differential equations (PDEs) using numerical methods. Some of the prevailing methods include the finite element method (FEM), mass-spring model (MSM), boundary element method (BEM), point-based method. Geometry-based models use intuitive methods instead of solving PDE. For example, vertex based model and spline based model. The advantages of geometry-based approaches are: fast and smooth, but the major drawback is in their accuracy. In contrast, physics-model is more accurate than geometry based approaches, but it has the disadvantages of being slow and rough. The key point to build a physical model of soft tissue is to reduce computational complexity and simultaneously keep the accuracy of deformation effects [2].

A material is called viscoelasticity if the stress–strain relationship of the material depends on time under applied forces. Ratchada [3] applied the method of tissue creep and stress relaxation to parallel-fibered collagenous tissues. Cotin [4] established a hybrid elastic model for surgery training and simulation. Suzuki [5] created a deformable organ model by using sphere-filled method instead of the finite element method. But the deformation of the model is still based on MSM. Recently, meshless methods have been introduced in the computer graphics field for a real-time simulation of deformable objects. The meshless method, in contrast to the mesh-based method, discretized the object by merely using a set of nodes without a mesh structure. This method thus avoided the problems associated with remeshing and handles the topology change of the object in a relatively simple manner. The topological relation among meshless simulation nodes is defined by support of the node. Adjacent nodes, which lie within the support, acted like the nodes which are connected by the edges of the mesh in the mesh-based method [6]. Then researchers start to set up the meshless physical model such as potential energy method [7], discrete mechanics [8], which has become a popular research in soft tissue modeling and cutting.

Haptic feedback is especially important in delicate surgical pressures requiring fine sensations, especially when slightly excessive pressure can injure a patient. The ability to touch the virtual objects in the environments is one factor that provides a more realistic experience and is particularly suited to virtual training applications [9]. So, the create process is related to the constraint force needed to keep a valid model which is close to the real tissue from the relaxation and creep and other properties.

MODEL STRUCTURE

Layered Sphere-Filled Structure

A multilayer structure has been adopted to arrange the TPM, as shown in Fig. (, ). The TPM has 424 spheres and the radius of the sphere is 10mm. These spheres are connected by a three-parameter structure which is composed of a Kelvin model and a linear spring.

Under external force, the deformation of the outer-layer sphere is obvious, and the deformation of the inner-layer sphere will weaken orderly. The deformation effects fit the deformation of the real soft tissue, which is shown in Fig. (, ). Fig. (, ) presents the structure of layered sphere-filled based on meshless method in TPM. This model has been established by using VC++ and the OpenGL.

To clarify this point further, the liver tissue model has been established by using MSM structure in this paper. Loading 2 newtons (N) in the Z-axis positive direction which act on the MSM, the deformation effects of the liver model are shown in Fig. (, ). Similarly, loading 2N in the Z-axis positive direction which act on TPM, the deformation effects of the liver model are shown in Fig. (, ). In the MSM model, the deformation effects of the force point and edges are roughness, as shown in Fig. (, , ). But TPM model has better effects in deformation, especially the deformation of the force point and edges.

Real-Time Simulation Condition

We used the hardware platforms that included Omega.7 force-feedback devices which came from the Swiss Force Dimension Company, a computer which consisted of AMD Athlon 64X2 2.0GHz CPU 2.0GB RAM and Samsung LCD, as shown in Fig. (, ). The graphic interfaces OpenGL can be used to display the effect of graphics rendering and visualization on this model simulation program. Programs of virtual scene setting and force-feedback effect came from the open source programs of the Swiss Force Dimension Company (http://www.chai3d.org /download.html).

MODEL FORMULATION

This model consisted of tension and three-parameter viscoelastic structure [10] which had creep effects and viscoelastic properties. Finally, stress relaxation and creep formulas had been obtained by mathematical formula derivation.

Model Framework

A three-parameter viscoelastic structure consists of a Kelvin model and a spring as the basic model was used in this paper, as shown in Fig. (). In Fig. (), we assume that is the stress of three-parameter structure of two nodes at t time, E and E are the elastic modulus of the spring, η is the viscous modulus, are the length of the spring at time.

It will be damaged when the deformation length of the liver is more than 12mm [11] by using a rigid object to touch the liver surface with different velocity. So, the deformation length of TPM model should remain within 10mm as an effective data by using a rigid sphere with a radius of 2mm to touch the model. Here, we set the radius of the sphere in TPM model as 10mm when exert external stress is exerted to the outer sphere in this model, the stress status as shown in Fig. (). It can be seen from Fig. () that the displacement of the outer sphere B changed under stress. To ensure the validity of the deformation, the displacement of sphere B should remain within 10mm. In order to simplify the calculation, we ignore the displacement of the spheres A and C. Here, we assume that is the center of the sphere A. Similarly, are the centers of spheres B, C, D, respectively. And is the center of the sphere B after deformation between two spheres which are connected by a three-parameter viscoelastic structure.

We assume that is the distance between and . The initial length of line AB is x, and the length of line AB is y after deformation. The angle between x and y is. Then, we assume . The length of deformation between A and B can be defined as

The deformation length of this model should remain within 12mm to avoid the model damage. In order to ensure this criterion, θ should satisfy that . Here, we assume that is the tension. Then, the tension can be defined as

where E' is the elastic modular.We assume , then Eq.(2) can be written as .

The formula of this model includes:

By using inverse Laplace transform, we can obtain

This assumes , , , Eq.(8) can be written as

Eq.(9) is the constitutive equation of this model.

Relaxation

After noting that . Here, is the unit step function. It can be defined as

By using Laplace transform , Eq.(9) can be written as

By using inverse Laplace transform, Eq.(10) can be written as

Eq.(11) is the relaxation equation of this model.

Creep

To exert a force on soft tissue, there is a step-response function . Considering load , substituting and into Eq.(9), by using Laplace transform, Eq.(9) can be written as

Eq.(12) can be written as

This assumes , by using inverse Laplace transform, Eq.(13) can be written as

Eq.(14) is the creep equation of this model.

Eq.(14) shows that the sphere will be stable under the constant force after t time, then the displacement of this sphere is . The relationship between and is the stress-strain in this model. The displacement ε is changed with time t that can be defined as creep. Eq.(11) shows that the force will change with the time when we give it a constant displacement . The force is changed with time t that can be defined as relaxation.

RESULTS

Model Parameters

The effective elastic modulus of soft organ was estimated from the force versus displacement data of static indentations based on the small deformation assumption. Carter [12] and Kauer [13] used an instrument manually to measure the displacement and force response of the organ surface. In experiments with ex vivo pig liver, by using a hand-held probe Carter [12] calculated that the elastic modulus of pig liver which was equal to 490 kPa . Ottensmeyer [14] designed a robotic indenter to measure mechanical properties of pig liver during a minimally invasive surgery. The robotic indenter can make small indentations in the range of ±500 µm. In experiments with in vivo pig liver, he calculated that the elastic modulus of pig liver was equal to 2.2–15 kPa.

Effective shear and elastic moduli of pig liver were estimated from the static indentation data by using the linear elastic contact theory and the small deformation assumption. The effective shear modulus, G, for small indentation of an elastic half space by a rigid hemispherical indenter can be obtained using the solution given by Lee [16]. So, by using the method of Lee the effective elastic modulus value was calculated between 11.6-18.7kPa [16]. On the other hand, in order to achieve better results, compared with the hand-held instrument, a robotic indenter was developed for minimally invasive measurement of the tissue properties. So we selected the parameter that was published by Ottensmeyer. Then we used the elastic modulus value of the liver which was equal to 12 kPa in this paper by combining these parameters.

It can be seen from Fig. (), under static stress loading, we have and . Then, we can write that .Therefore, we get that =12kPa. For ,, the value of is calculated between 0.096kPa and 5.9kPa.

Klatt [17] had concluded that the best agreement between fit and experiment was achieved with the Zener model. The reported value for viscous modulus base on Zener viscoelastic model was equal to 5.5 ± 1.6 Pa.s. Mohammad [18] calculated that the viscous modulus value was between 0.8-1.57 Pa.s. Asbach [19] calculated that the viscous modulus value was equal to 5Pa.s. This deviation between our modeling results and other experiment results not only because of the model type (Voigt or Zener), but also because of different experimental conditions such as the slice thickness, kind of liver and in vivo or ex vivo execution. In this paper, the structure of TPM model was similar to that reported in Asbach [19], so we chose to use 5Pa.s as our viscous modulus value in TPM.

Model simulation

A force-feedback instrument was used to touch the No.71(x=0.002, y=0.136, z=0.057) sphere in the simulation experiments. The motion of other sphere was driven by the variation of the No.71 sphere. And we selected above-mentioned parameters as the appropriate parameters for constitutive equations.

Evren [15] used a probe with a diameter of 4mm to insert into the abdominal cavity of three different pigs. In the static loading test, the liver of each pig will indent slowly from 0 to 8mm to characterize its force and displacement response. In Fig. (), the blue curve, green curve and cyan curve represent the force response value of the liver of pig #1, pig #2 and pig #3, respectively. It can be seen from Fig. () that the force response value of three pigs for the indentation depth is different, but these data include a valid data range. In order to match the model experiment, a sphere with a diameter of 4mm was used as a force-feedback instrument in TPM model, which is shown in Fig. (, , ). As shown in the Fig. (), red curve represents the force response value of No.71 sphere in TPM for the indentation depth changing from 0 to 8mm. The force response value of TPM model for the indentation depth maintains between the pig #1 and the pig #3.

Noting that when in Eq. (11), we have

Here, we set that .

Similarly, noting that when, we have

If the indentation depth of this model is 2 mm, as shown in the Fig. (), the distance of two spheres is x, which will influence the model precision. When we assume that x is equal to 20mm, =0.119kPa can be obtained by using Eq. (2). Substituting into the Eq (11)，the relaxation data in TPM model can be obtained. Similarly, when x is 25mm, 30mm and 35mm respectively, will change as the value of x. It can be seen from Table , the result shows that both the model simulation data and the experimental data are very close. In the same way, when the indentation depth of this model is 4mm, 6mm and 8mm respectively, the relaxation data in TPM model is compared with the real liver experimental data, as shown in the Table , Table and Table .

By analyzing these tables, we can find that the relaxation data in TPM are close to experimental results when x=30mm. Because the radius of the sphere in TPM is 10mm, so we set that the distance of two spheres is 10mm.

Here, we assume x=30mm and =4mm, then the angle between two outer spheres should satisfy , and we can get . Substituting θ into the Eq.(2) =0.211kPa can be obtained. Finally, using , we can obtain F=0.613N, which means the instantaneous forces of liver. The relaxation curve is recorded at the indentation depth of 4 mm at strain rates of 4mm/s, 2mm/s, and 1 mm/s respectively in TPM, as shown in Fig. ().

Table shows the comparison between the force response value of TPM and a valid data range of the force response value which was obtained from pig liver experiments. It can be seen from Tab.5 that the force response value of TPM is in the range of the pig liver experimental data. So, it can be verified that the force response value in TPM was effective.

In order to further verify the effectiveness of TPM, we compared the relaxation data in TPM with the relaxation data that was reported in Evren [15]. In stress relaxation experiments, the liver of each pig was indented to pre-defined depths of 2mm, 4mm, 5mm, 6mm, and 8mm respectively and recorded the force relaxation response value in the next 30 seconds. In the Fig. (), the blue curve, green curve and cyan curve represent the stress relaxation response value of pig #1, pig #2 and pig #3, respectively. The red curve represents the stress relaxation response value of the #71 sphere in TPM at the indentation depth of 4mm with a strain rate of 4mm/s in the next 30 seconds. It can be seen from Fig. () that the lines of stress relaxation value of TPM are close to the experimental data of pig liver.

Here, we select that the force acting on the model is 0.196N. We assume λ=0.179kPa. Noting that when , we have

Then we can obtain the indentation depth of the model which is 1.31mm by substituting these parameters into Eq.(17). Similarly, noting that when , we have

Then we can obtain the indentation depth of this model which is 3.8mm by substituting these parameters into Eq.(18). Similarly, when we select that the force acting on the model is 0.96N and we assume λ=5.8kPa by using Eq.(6). Noting that when , we can obtain the indentation depth of TPM which is 4.4mm by substituting these parameters into Eq.(17). Noting that when , we can obtain the indentation depth of TPM which is 10.86mm by substituting these parameters into Eq.(18).

In order to further verify the effectiveness of TPM, we compared the creep data in TPM with the creep data that was reported in Amy [20] .In the Fig. (), the cyan curve and the blue curve represent the lower boundary and upper boundary of creep data of the Ex vivo perfused pig liver under 0.196N respectively. The indentation depth of the Ex vivo perfused pig liver maintains between 3.8mm and 6.65mm. The red curve represents the creep data of the #71 sphere in TPM under 0.196N. It can be seen from the Fig. () that the creep curve in TPM is in the range of lower boundary to upper boundary. In the Fig. (), the cyan curve and the blue curve represent the lower boundary and upper boundary of creep data of the Ex vivo perfused pig liver under 0.98N respectively. The indentation depth of the Ex vivo perfused pig liver maintains between 9.6mm and 11.9mm. The red curve represents the creep data of the #71 sphere in TPM under 0.98N. It can be seen from Fig. () that the creep curve in TPM is in the range of lower boundary to upper boundary.

DISCUSSION AND CONCLUSIONS

Comparing with the experimental data of pig liver which was reported in Evren [15] and Amy [20], we verified that the TPM was effective in stress-strain, stress relaxation and creep data. Moreover, the TPM with two-layer spheres structure had been established. On the basis of this model, in order to further improve the model precision, we can establish TPM structure with three or several layer spheres.

Although this model had stress-strain, stress relaxation and creep deformation properties, which had a significant advantage over the MSM, but the precision of this model needs to be improved further. And the data gained from the simulation experiment had some deviation. In the future, we will refine the model to increase the precision of the model, and the model will be evaluated on more challenging and complicated tasks.

Table 1.

When the Indentation Depth is 2 mm, the Different Values of λ Were Obtained with Different x. The Table Shows the Comparison between The Relaxation Data of TPM Model and the Real Liver Experimental

x (mm)	(kPa)	The Relaxation Data of TPM Within 30s (N)	The Real Liver Experimental Relaxation Data Within 30s(N)
20	0.119	0.304→0.153	0.32→0.21
25	0.076	0.303→0.152	0.32→0.21
30	0.053	0.302→0.152	0.32→0.21
35	0.039	0.301→0.151	0.32→0.21

Table 2.

When the Indentation Depth is 4 mm, the Different Values of λ were Obtained with Different x. The Table Shows the Comparison between the Relaxation Data of TPM Model and the Real Liver Experimental

x(mm)	(kPa)	The Relaxation Data of TPM Within 30s (N)	The Real Liver Experimental Relaxation Data Within 30s(N)
20	0.466	0.626→0.324	0.62→0.32
25	0.301	0.618→0.316	0.62→0.32
30	0.211	0.613→0.312	0.62→0.32
35	0.155	0.61→0.3	0.62→0.32

Table 3.

When the Indentation Depth is 6 mm, the Different Values of λ were Obtained with Different x. The Table Shows the Comparison between the Relaxation Data of TPM Model and the Real Liver Experimental

x(mm)	(kPa)	The Relaxation Data of TPM Within 30s (N)	The Real Liver Experimental Relaxation Data Within 30s(N)
20	1.01	0.98→0.528	0.88→0.48
25	0.662	0.95→0.5	0.88→0.48
30	0.466	0.93→0.48	0.88→0.48
35	0.345	0.928→0.47	0.88→0.48

Table 4.

When the Indentation Depth is 8 mm, the Different Values of λ were Obtained with Different x. The Table Shows the Comparison between the Relaxation Data of TPM Model and the Real Liver Experimental

x(mm)	(kPa)	The Relaxation Data of TPM Within 30s (N)	The Real Liver Experimental Relaxation Data Within 30s(N)
20	1.71	1.37→0.774	1.72→0.68
25	1.14	1.32→0.71	1.72→0.68
30	0.81	1.28→0.68	1.72→0.68
35	0.6	1.26→0.66	1.72→0.68

Table 5.

The Comparison Between the Force Response of TPM Model and a Valid Data Range of the Force Response was Obtained Through Pig Liver Experimental

(mm)	(kPa)	The Value of Force in TPM Model (N)	A Valid Data Range of the Force (N)
2	0.053	0.152	0.07～0.18
4	0.211	0.312	0.13～0.36
6	0.466	0.48	0.2～0.6
8	0.81	0.68	0.23～0.89